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ABSTRACT 


The  basic  filter-observer  equations  of  Kalman  for  optimal  and 
suboptimal  filters  are  studied  usinq  the  concepts  of  Lyapunov  functions 
and  stability  theory.  The  Second  Method  of  Lyapunov  is  used  to  form  a 
basis  for  comparison  of  the  convergence  rates  of  such  filters.  Lyapunov 
functions  are  also  used  to  derive  constraining  relations  for  the 
elements  of  the  filter  gain  matrix  leading  to  design  criteria  for  sub- 
optimal  filters.  A  derivation  of  the  optimal  filter  gain  based  upon  the 
Lyapunov  function  of  a  random  variable  is  given.  This  derivation  shows 
that  the  optimal  filter  converges  most  rapidly.  A  desiqn  of  a  sub- 
optimal  filter  for  one  class  of  signal  models  is  given  based  solely  upon 
stability  constraints. 
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SYMBOLS  AND  ABBREVIATIONS 

X         an  n  square  matrix  unless  otherwise  specified 

x         an  n  vector.  Scalars  will  be  specified  unless  obvious 
from  the  notation 

t         independent  variable  time  and  is  always  a  scalar 
(similarily  for  n  in  discrete  time) 

det(A)  determinant  of  matrix  A 

tr(A)  trace  of  matrix  A 

A (A)  eigenvalue  of  matrix  A 

|-||  norm  of  the  argument  (see  Appendix  A) 

| ot |  magnitude  of  the  complex  scalar  a 

N(x  ,6)     neighborhood  about  the  point  x   of  radius  6  . 

Taken  as  the  usual  sphere  definition  all  x  such  that 

I|X  -  Xfi||   1  6 

E[x]       expected  value  of  the  random  vector  x 
e         is  a  member  of,  belongs  to 
inf[«]     greatest  lower  bound 
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I.   INTRODUCTION 

For  many  problems  in  control  it  is  desirable  to  have  some  measure 
of  all  states  in  a  given  signal  process.  By  their  nature,  many 

processes  only  allow  us  to  physically  measure  or  observe  some  of  the 
states  characterizing  the  process.  Moreover,  there  is  some  amount  of 
uncertainty  in  the  measured  states  which  may  be  characterized  by 
additive  measurement  noise.  Thus,  estimation  theory,  as  a  device  for 
estimating  the  states  of  a  system  in  the  presence  of  noisy  and  incom- 
plete observations,  has  emerged  as  an  important  part  of  modern  control 
theory. 

In  the  early  1960 ' s  Kalman  [K3,K4]  advanced  the  well  known  Wiener 
filter  theory  [ref.  D2,P1]  by  showing  how  such  linear  least  square 
filters,  which  Wiener  synthesized  using  classical  frequency  domain 
theory,  may  be  realized  in  the  time  domain.  This  development  coin- 
cided with  the  development  of  the  digital  computer  and  its  use  as  an 
active  component  in  process  control  systems. 

The  basic  filter  equation  of  Kalman  is  really  a  recursive  weighted 
average  technique  of  the  following  form: 

X  =  X  +  G  (X-Z)  (1) 

where  X  is  the  corrected  estimate,  X  is  the  projected  estimate 
following  known  signal  dynamics,  G  is  the  filter  gain  weighting 
factor,  and  Z  is  the  signal  observation  including  noise.  G  is 
generally  time  varying  and  its  selection  constitutes  the  basic  design 
and  nature  of  the  filter.  The  difference  equation  of  (1)  becomes  a 
differential  equation  in  continuous  time. 


The  corrected  estimate,  X  ,  is  usually  an  n  vector  of  the  states 
of  the  signal  process,  and  the  observation,  Z  ,  is  an  m  vector  with 
m  <_  n. 

The  selection  of  the  filter  gain  has  two  objectives: 

1.  To  provide  an  estimate  of  unobserved  states 
based  upon  data  from  the  observed  states. 

2.  To  provide  for  some  weighting  of  the  observations 
to  allow  compensation  for  the  measurement  noise. 

Objective  1  is  the  basic  observer  problem  which  has  emerged  in  recent 

years  (such  as  the  Luenberger  observer  [LI])  and  does  not  consider 

any  effects  of  noisy  measurements.  Objective  2,  as  stated,  involves 

estimation  theory  and  the  stochastic  properties  of  the  observation 

noise  as  well  as  the  random  input  signals  which  drive  the  system. 

Current  literature  contains  many  examples  of  the  application  of 
Kalman  filters  to  various  problems.  The  computational  capability 
required  to  implement  such  filters  is  generally  large,  and  often  a 
much  simpler  implementation  works  nearly  as  well.  The  literature  also 
indicates  that  the  method  generally  employed  to  determine  the  parame- 
ters of  these  filter  implementations  involves  large  scale  simulation 
trials. 

This  thesis  studies  the  properties  of  the  basic  filter  observer 
equation  in  both  continuous  and  discrete  time  using  the  Second  Method 
of  Lyapunov.  Lyapunov  theory  [ref.  K1,H1]  was  developed  for  investi- 
gating the  stability  properties  of  systems.  Its  application  to  the 
filter-observer  problem  has  appeared  in  recent  literature  (Deyst  and 
Price  [Dl]).  However,  this  approach  has  not  been  fully  exploited  in 
the  literature,  and  its  possibilities  are  investigated  in  this  thesis. 
This  study  leads  to  engineering  insight  into  the  operation  of  suboptimal 
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filters  and  provides  some  basis  for  comparison  of  the  performance  of 
these  filters  when  applied  to  the  same  problem.  Some  specific  results 
obtained  include  constraining  relations  on  the  filter  gain  matrix 
which  insure: 

1.  that  the  filter  is  asymptotically  stable. 

2.  that  the  filter  will  converge  to  steady  state 
at  a  rate  faster  than  a  given  bound. 

Application  of  the  filter  observer  equations  as  discussed  here  is 
limited  to  linear  systems  and  it  is  assumed  that  the  signal  dynamics 
are  known  exactly.  A  review  of  stability  concepts  for  linear  time- 
varying  systems  and  the  Second  Method  of  Lyapunov  is  included. 

Chapter  II  reviews  basic  definitions  and  stability  theorems  for 
linear  systems,  based  upon  the  concepts  of  state-space  theory,  with 
emphasis  on  forced  systems  so  that  these  results  may  be  applied  readily 
to  the  filter-observer  problem. 

Chapter  III  reviews  the  Second  Method  of  Lyapunov.  This  chapter 
includes  some  theorems  on  the  determination  of  Lyapunov  functions  for 
certain  classes  of  systems  as  well  as  some  theorems  on  transient  esti- 
mation using  Lyapunov  theory  which  have  not  appeared  in  the  literature 
previously. 

Chapter  IV  applies  stability  theory  to  the  filter-observer  homo- 
geneous dynamics  and  contains  many  original  results  including  the 
formulation  of  a  time  invarient  Lyapunov  function  for  the  filter.  The 
same  Lyapunov  function  may  be  used  for  a  large  class  of  filters 
applied  to  the  same  problem.  Particular  filters  have  different  time 
derivatives  of  the  Lyapunov  function.  This  leads  to  a  method  of 
comparing  convergence  rates  for  various  filters  applied  to  the  same 
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problem.  By  applying  the  Hadamard-Gerschgorin  theorem  to  provide 
sufficient  conditions  for  positive  definiteness  of  a  symmetric  matrix, 
certain  constraining  relations  are  derived  for  the  elements  of  the 
filter  gain  matrix  which  lead  to  the  design  of  stable  filters. 

Chapter  V  considers  the  forced  filter  dynamics.  New  results 
include  the  introduction  of  the  concept  of  a  Lyapunov  function  for  a 
random  variable.  A  derivation,  using  these  Lyapunov  functions,  of 
the  optimal  filter  is  given  showing  that  such  filters  also  converge 
most  rapidly  to  the  minimum  covariance  of  error. 

Chapter  VI  applies  the  results  of  Chapter  IV  to  the  design  of  a 

J-  L. 

class  of  suboptimal  filters  for  the  general  n   order  signal  model 
in  phase  variable  form. 
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II.  STABILITY  OF-  LINEAR  SYSTEMS 

A.  INTRODUCTION 

A  review  of  stability  definitions  and  theorems  are  set  forth  in  the 

following*  Many  current  references  are  available  for  this  material 

(HI,  Kl ,  K2,  01,  SI,  Tl  and  others).  Although  here  we  are  mainly 

interested  in  linear  systems,  many  of  the  theorems,  as  noted,  apply 

to  nonlinear  systems .  The  reader  is  referred  to  Appendix  A  for  a 

summary  of  norms  and  norm  properties. 

Definition  2.1    Free  system:  Any  system  with  no  input 

forcing  functions.  Such  a  system  is 
described  as  x  =  f(x,t) 

Definition  2.2    Autonomous  System:  A  free  system  which 

is  also  time  invariant,  (i.e.,  x  =  f(x)) 

Definition  2.3    Equilibrium  state (x  ):  Any  state  for 


which  x  =  0.  For  a  free  system 
f(xe,t)  =  0 

Definition  2.4    Boundedness:  An  equilibrium  state  x 


e 
is  said  to  be  bounded  if,  and  only  if, 

there  exists  a  neighborhood  N(x  ,6) 

such  that,  if  the  initial  state  x(t  ) 

lies  within  N,  then  ||x(t)-x  ||  <_  m  < 

for  all  t  >  t 
o 

Definition  2.5    Stability  (in  the  sense  of  Lyapunov): 
An  equilibrium  state  x  is  said  to  be 
stable,  if,  and  only  if,  to  any  neigh- 
borhood N(x  ,e)  there  corresponds  a 
neighborhood  N(x  ,6)  such  that,  if 
x(t  )  lies  in  N(x  ,6)  then  x(t)  lies 

in  N(x  ,t)  for  all  t  >  t  . 
e  o 


13 


The  difference  between  the  last  two  definitions  is  important  for 

nonlinear  systems  although  the  definitions  are  equivalent  for  linear 

systems.  Note  that  Definition  2.5  requires  that  the  state  remain 

within  a  preassigned  neighborhood  (which  may  be  arbitrarily  small)  for 

all  t  >  t  whereas  the  definition  of  boundedness  requires  only  that  the 

distance  (in  state  space)  from  x  to  the  trajectory  x(t)  remain  bounded, 

Thus,  boundedness  allows  for  such  things  as  limit  cycles. 

Definition  2.6    Asymptotic  Stability:  An  equilibrium  state 

is  said  to  be  asymptotically  stable  if 

(1)  It  is  stable. 

(2)  |  |x(t)-x  |  |  approaches  zero  as  t 
approaches  infinity  for  any  initial  con- 
dition lying  in  the  neighborhood  N(x  , 6) 
for  which  it  is  stable. 

The  above  definitions  are  for  strictly  local  conditions  about  the 
equilibrium  state.  If  the  neighborhood  N(x  , 6)  can  be  the  entire  state 
space,  then  the  system  is  said  to  be  globally  stable  or  asymptotically 
stable  in  the  large  (ASIL)  about  the  equilibrium  state  x  . 

Uniform  stability  is  also  a  useful  concept.  It  means  that  the 
neighborhood  in  which  the  initial  state  must  lie  is  independent  of  the 
initial  time.  That  is,  <5  is  not  a  function  of  t  . 

Definitions  2.1  through  2.6  above  alsu  apply  to  discrete  time 
systems  with  the  use  of  (n,n  )  in  place  of  (t,t  )  respectively.  In 
the  sequel,  theorems  and  concepts  are  developed  for  both  continuous 
and  discrete  time  systems.  The  discrete  notation  is  as  follows: 


Some  authors  distinguish  between  global  and  stability  in  the 
large  (HI)  but  for  linear  systems,  they  are  identical. 
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The  System: 

x(n+l)  =  AD(n)  x(n)  +  BD(n)  u(n)  (2.1) 

The  Observation: 

y(n)  =  HD(n)  x(n)  (2.2) 

The  Fundamental  Matrix: 

*(n,nQ)  =  AD(n-l)  AQ(n-2)  .  .  .  AD(nQ)  (2.3) 

The  Solution: 

n-1 

x(n)  ■  *(n,n  )  x(nj  +  z     *(n,k+l)  B(k)  u(k)      (2.4) 

0    °   k=n 
o 

For  time  invariant  systems  we  have  $(n)  =  ADn  . 

Also  it  is  interesting  to  note  the  relationship  between  continuous 

sampled  systems  and  their  discrete  equivalent, for  the  time  invariant 

case. 

T 
AD  =  $(T)     BD  =   /  $(T-t)  B  dx  (2.5) 

o 

B.   FREE  SYSTEMS 

It  is  well  known  that  for  linear  time  invariant  systems,  if  the 
eigenvalues  of  the  system  matrix  A  all  have  negative  real  parts,  then 
the  system  is  asymptotically  stable.  Similarly  for  discrete  systems 
the  system  is  asymptotically  stable  to  the  origin  if,  and  only  if,  the 
magnitude  of  the  eigenvalues  of  AD  are  less  than  unity.  Such  condi- 
tions are  not  easily  established  for  time  varying  systems.  However, 
the  following  theorems  apply  for  the  system  x(t)  =  A(t)  x(t),  with 
x  =  0  and  the  equivalent  discrete  case  [ref.  SI]. 

Theorem  2.1      The  origin  of  a  linear  system  is  bounded 

if,  and  only  if,  the  norm  of  the  fundamental 

matrix  <J>(t,t  )  (or  $(n,n  ))  is  bounded. 
o         o 
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To  prove  sufficiency  note  that 

||x(t)||  =  ||*(t,t0)  x(tQ)  It    £   Il*(t,t0)||  ||x(to)|| 

If      ll*(t,t0)||  <  K 

Then    ||x(t)||  <  K||x(t0)|| 

By  definition  2.4,  with  6  =  ||x(t  )||   it  follows  that  the  origin 
(x  =0)  is  bounded.  Necessity  is  shown  by  considering 
||x(t)||  =  ||*(t,t0)  x(to)|| 

and  noting  that  if  the  norm  of  *(t,t  )  is  not  bounded  then  some 

element  of  *(t,t  )  approaches  infinity  with  time  and  hence  ||x(t)| 

also  approaches  infinity  and  is  unbounded. 

Theorem  2.2      Boundedness  and  stability  of  the 

origin  are  equivalent  in  linear  systems. 

To  prove  this,  we  must  show  that  for  any  e  >  0  it  is  possible  to  find 

a  6  such  that  if  x(t  )  is  in  N(0,6)  then  x(t)  is  in  N(0,e)  for  all 

t  >  t  .  This  reduces  to  showing  that  if  ||x(t  )||  <  6  then  ||x(t)||  <  e 

By  definition  2.4,  if  the  origin  is  bounded  we  may  take  6  =  -|—  . 

Then  by  theorem  2.1 

Mx(t)||  1K||x(to)||  <  J%-  =  e 


which  shows  that  boundedness  implies  stability.  It  is  obvious  that 

stability  implies  boundedness  hence  the  two  are  equivalent. 

Theorem  2.3      The  origin  of  a  linear  system  is 

asymptotically  stable  if,  and  only  if, 
the  norm  of  $(t,t  )  is  bounded  for  all 
t  <  t  and  $(t,t  )  approaches  zero  as 
t  approaches  infinity. 
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If  ||*(t,t  ||  ^K  then  the  origin  is  stable  (Theorems  2.1  and  2.2).  We 


o 
now  have  to  show  that  ||x(t)|  approaches  zero  as  t  approaches  infinity. 

Now  ||x(t)||  =  ||*(t,t  )  x(t  )||  and  if  *(t,t  )  approaches  zero  then 

||x(t)||  approaches  zero. 

Moreover,  if  ||x(t)||  approaches  zero  for  all  ||x(t  )||  in  N(0,6) 
then  $(t,t  )  approaches  zero,  since  x(t  )  is  arbitrary  in  N(0,6). 

It  should  be  noted  that  if  a  linear  system  is  asymptotically  stable, 
then  it  is  also  ASIL,  since  the  conditions  for  stability  depend  only 
upon  the  fundamental  matrix  and  are  independent  of  the  initial  condi- 
tions of  the  states. 

The  foregoing  theorems  also  apply  to  discrete  systems  with  the 
change  in  notation  from  <J>(t,t  )  to  $(n,n  ). 

C.   FORCED  SYSTEMS 

Now  we  are  interested  in  studying  forced  systems  of  the  form 

x(t)  -  A(t)  x(t)  +  B(t)  u(t)  (2.6) 

We  are  particularly  interested  in  the  effect  of  u(t)  on  the  system, 
assuming  it  is  at  rest  at  t=t  .  Then  the  solution  of  (2.6)  is 

x(t)  =   /  *(T,t_)  B(t)  u(t)  dT  (2.7) 

Stability  of  forced  systems  is  considered  with  respect  to  a  given 

set  of  inputs  U  =  {u(t)}  .  Generally,  this  set  is  taken  to  be  the 

set  of  bounded  inputs. 

Definition  2.7    Stability  of  a  forced  system. 

A  forced  system  is  stable  with 
respect  to  U  if,  and  only  if, 
the  states  are  bounded  (||x(t)|| 
<  K)  for  all  u(t)  in  U  and  all 
t  >  t  . 
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Theorem  2.4      A  forced  system  is  stable  with 

respect  to  the  set  of  bounded  inputs 

ift 

f.    ll*(T,t0)  B(T)||  di  <  K1 
o 

To  prove  sufficiency  note  that 

t  t 

l|x(t)||  <  /  ||*(T,t. )B(t)u(t)||  di  <  /  ||»(T,t  )B(x) | |  ||u(T)||dT 

tn  °  t        ° 

o  0 

If  I  |u(t)||  <  K2    then  |  |x(t)|  |  <  K^  . 

Hence  the  state  is  bounded. 

A  consolidating  theorem  appears  in  Timothy  and  Bona  [Tl]  and 

Kalman  [Kl]  with  proof  given  in  [Kl]. 

Theorem  2.5      For  the  system  described  by  equation 

(2.6)  if 

(i)  | |A(t)| |  <  K1  for  all  t 

(ii)  0  <  K2  <  ||B(t)v||  <  K3  for  all   ||v||=l 

where  v  is  an  arbitrary  vector, 
then  the  following  statements  are 
equivalent. 

1.  Any  uniformly  bounded  input 
I lu(t) I  I  1  K4  results  in  a 
bounded  output  ||x(t)||  <  Kr  <  °° 

2.  For  all  t  >  t 

0 

t 

/   ||*(t,T)||dxl  K6  <  - 

to 

3.  The  free  system  is  uniformly 
asymptotically  stable. 
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III.  THE  SECOND  METHOD  OF  LYAPUNOV 

A.   INTRODUCTION 

The  foregoing  indicates  that  stability  for  linear  systems  is 

completely  determined  by  the  boundedness  properties  of  the  fundamental 

matrix  and  the  inputs  to  the  system.  This  means  in  essence  that  the 

complete  solution  to  the  system  must  be  obtained  before  stability  can 

be  ascertained.  Lyapunov  proposed  a  method  by  which  stability  of  a 

system  may  be  determined  without  finding  the  solution  to  the  system 

equations.  His  method  is  based  upon  a  generalized  energy  concept. 

From  physics  it  is  known  that  for  a  nonconservative  system,  if  the 

total  energy  is  always  decreasing,  then  the  system  is  stable.  The 

total  energy  of  such  a  system  is  a  scalar  function  of  the  state  of  the 

system.  Lyapunov  has  shown  that  the  very  existence  of  scalar  functions 

of  the  state  of  a  system  which  obey  certain  properties  is  sufficient 

to  conclude  that  the  system  is  stable.  These  functions  are  called 

Lyapunov  functions.  It  is  interesting  to  note  that  Lyapunov  functions 

for  a  stable  system  are  not  unique  and  that  a  sufficient  condition  for 

stability  is  the  existence  of  any  such  function. 

Definition  3.1    Lyapunov  function.  Any  function,  V(x), 

having  the  following  properties  is 
called  a  Lyapunov  function. 

(a)  V(x)  is  a  continuous  scalar  function 
of  the  systems  state  vector  and  has 
continuous  first  partial  derivatives. 

(b)  V(x)  is  positive  definite. 

(c)  V(x)  =  [grad  V(x)]  x  is  negative 
definite. 
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Definition  3.2    A  scalar  function,  V(x),  of  a  vector  argu- 
ment is  positive  definite  if, and  only  if, 

(a)  V(0)  =  0 

(b)  V(x)  >  0   whenever  x  f   0 

Positive  semi-definite  means  that  equality 
is  also  allowed  in  part  (b). 

For  negative  definite  functions  the 
inequality  is  reversed. 

In  two  dimensional  state  space  the  Lyapunov  function  is  a  cupped 

shaped  surface  setting  on  the  origin.  Its  continuity  and  definiteness 

mean  that  the  level  curves  of  V  projected  on  the  state  plane  are 

closed  about  the  origin  and  the  curve  k,  =  V(x)  is  inclosed  within 

k2  =  V(x)  if  k,  <  k2  . 

Under  the  foregoing  conditions  we  can  consider  that  the  Lyapunov 

2 
function  gives  a  measure  of  distance  from  the  origin  in  state  space 

as  the  system  follows  the  trajectory  x(t). 

B.  THE  BASIC  STABILITY  THEOREM 

Before  stating  the  theorem,  consider  the  effect  of  explicit  time 
variation  upon  the  properties  of  the  Lyapunov  function.  Suppose  we 
have  a  function  V(x,t)  and  its  time  derivative 

V(x,t)  =  -U-  +  [grad  V]T  x  (3.1) 

Now,  V(x,t)  will  be  positive  for  all  t  if  it  is  always  greater 
than  a  positive  function  e(||x||).  However,  this  is  not  sufficient 
to  guarantee  asymptotic  stability  even  if  9  <  0  for  all  t  and  x  f   0. 


2 
Distance  used  in  this  sense  is  not  the  same  as  the  usual 

Euclidian  distance.  Distance  is  implied  by  the  ordered  nature  of  the 

level  curves  of  V  as  discussed. 
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Since  the  grad  V  term  in  equation  (3.1)  indicates  the  motion  of  the 

system  state  (i.e.,  is  multiplied  by  x),  it  is  possible  that  -|i- 

be  sufficiently  negative  so  that  v"  <  0  for  all  t,  while  the  state 

moves  outside  of  any  bounding  region  (N(x  ,e)).  Hence,  the  state  may 

never  go  to  the  equilibrium  state  x  .  Hsu  and  Meyer  [HI]  have  a  good 

discussion  and  Kalman  and  Bertram  [Kl]  provide  a  rigorous  proof  of 

the  following  theorem. 

Theorem  3.1      Consider  the  continuous  time,  free 

dynamics  system  x  =  f(x,t)  where  f(0,t)  =  0 
for  all  t.  Suppose  there  exists  a  scalar 
function  V(x,t)  with  continuous  first 
partial  derivatives  such  that  V(0,t)  =  0  and 
(i)  there  exists  a  continuous,  non- 
decreasing  scalar  function  a  such  that 
a(0)  =  0  and  V(x,t)  >_  ct(  |  |x|  | )  for  all  t 
and  x  f   0. 

(ii)  there  exists  a  continuous  positive 
scalar  function  y   such  that  y(0)  =  0  and 
V(x,t)  i-y(||x||) 

(iii)  v(x,t)  <  e(||x||) 

where  s  is  a  continuous  nondecreasing 
function  such  that  (3(0)  =  0  for  all  t. 

(iv)  a(  |  |x|  | )  -*  °°  with  |  |x|  |  +  °° 

Then  the  system  is  asymptotically  stable  in  the  large  to  the 

equilibrium  state  x  =  0  and  V(x,t)  is  a  Lyapunov  function  for  the 

system. 

Requirement  (i)  insures  that  the  Lyapunov  function  is  always 

positive  whereas  (iii)  insures  that  the  function  does  not  stay  infinitely 

large  throughout  the  state  space.  In  particular  it  is  required  that  the 

V  function  does, in  fact, go  to  zero  uniformly  with  time  at  the  origin. 
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The  y  function  in  (ii)  assures  that  V  is  always  negative.  This 
together  with  (1)  and  (iii)  assure  that  the  Lyapunov  function  following 
the  system  trajectory  is  decreasing  and  must  end  up  at  the  origin. 
To  insure  that  the  conditions  exist  for  all  of  state  space,  require- 
ment (iv)  is  imposed  and  the  system  is  ASIL.  It  should  be  noted  that 
a  positive  definite  quadratic  form  for  V  automatically  satisfies  the 
requirements  of  (i),  (iii),  and  (iv). 

For  discrete  systems  we  define  a  Lyapunov  function  as  V(x,n) 
and  the  rate  of  increase  of  V  along  the  system  trajectory  as 

AV(x,n)  =  V(x(n+l),n+l)  -  V(x(n),n)  (3.2). 

With  these  changes,  Theorem  3.1  becomes  the  corresponding  theorem 
for  discrete  systems. 

C.  THE  LINEAR  AUTONOMOUS  SYSTEM 

For  linear  autonomous  systems,  it  is  a  necessary  and  sufficient 

condition  for  ASIL  that  a  quadratic  form  Lyapunov  function  exist. 

This  is  embodied  in  a  theorem  originally  proved  by  Lyapunov. 

Theorem  3.2      The  origin  of  a  linear  autonomous 

system,  is  asymptotically  stable  if, 
and  only  if,  for  any  symmetric  positive- 
definite  matrix,  C,  there  exists  a 
unique  positive-definite  matrix  Q  which 
satisfies  the  matrix  equation 

ATQ  +  OA  =  -C  (3.3) 

A  full  proof  of  this  theorem  appears  in  Hahn  [H2].  Essential  parts 

also  are  shown  by  Bellman  [Bl].  The  major  results  of  interest  are 

delineated  here  for  convenience. 

1.  The  Lyapunov  function  is  V  =  x  Qx 
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2.  By  direct  computation  , 

tf  =  xTQx  +  xTQx  =  xTATQx+  xTQAX 

T  (3.4) 

=  -x'Cx 

3o  Given  any  positive  definite  matrix  C  then  Q  is  uniquely 
determined  by  the  matrix  equation  (3.3)  if  no  two  eigen- 
values of  A  sum  to  zero  or  no  eigenvalue  is  zero,  (see  [Bl].) 

4.  If  A  is  stable  then  the  unique  Q  is  given  by 

00  ATt    At 
Q  =   /  eM  z     C  eML  dt  (3.5) 

0 

3 
Using  4  and  a  theorem  stated  in  Browne  [B2]  it  is  easy  to  show 

that  if  C  is  positive  definite,  then  0  must  be  also.  Since  C  =  B  B 

3 
where  B  is  a  real  non-singular  matrix  we  can  write 

0  =   /  eA  tBTBeAt  dt  =   /  (BeAt)  (BeAt)  dt  (3.6) 

0  0 

At  At 

But  B  and  e   are  non-singular.  Hence  Be  is  non-singular  and  0  is 

3 
positive  definite  by  reappli cation  of  Browne's  theorem  . 

Generally  the  use  of  Theorem  3.2  proceeds  as  follows: 

(i)  Choose  any  matrix  C  which  is  positive  definite  and  solve  the 

n(n+l)/2  algebraic  equations  (3.3)  for  the  elements  of  Q. 

(ii)  Test  Q  for  definiteness.  If  it  is  positive  definite  the  system 

is  asymptotically  stable.  If  it  is  not  positive  definite  the  system 

is  not  asymptotically  stable.  This  follows  from  the  necessary  and 

sufficient  conditions  stated  in  Theorem  3.2. 


3  T 

From  Browne;  If  C  is  any  real  non-singular  matrix  and  A  =  C  C 

then  A  is  positive  definite.  The  converse  is  true,  any  positive 

definite  matrix  A  can  be  expressed  as  such  a  product.  Similarily,  if 

C  is  n  square  of  rank  r<n  then  A  is  positive  semi-definite  of  rank  r. 
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The  theorem  corresponding  to  Theorem  3.2  for  discrete  time  systems 

is  as  follows  [ref.  Kl,01]. 

Theorem  3. 2D     The  discrete  time  system  described  by 

x(n+l)  =  AD  x(n) 

is  asymptotically  stable  if,  and  only  if, 
for  any  given  positive  definite  matrix  C 
there  exists  a  unique  positive  definite 
matrix  Q  satisfying 
T, 


AD'QAD  -  Q  =  -C 
where  V  =  x  Qx  is  a  Lyapunov  function 


(3.7) 


with  AV  =  -x  Cx. 

An  algorithm  for  numerical  solution  of  equation  (3.7)  for  Q  given 
any  positive  definite  matrix  C  is  derived  in  Appendix  B.  This  is  a 
new  result  and  makes  the  application  of  Theorem  3. 2D  much  easier  than 
in  the  past. 

At  this  point  it  is  important  to  note  that  we  cannot  in  general 
choose  any  positive  definite  matrix  for  Q  and  then  perform  the  multipli 
cation  and  addition  indicated  in  (3.3)  to  obtain  a  C  matrix  which  is 
positive  definite,  even  if  A  is  known  to  be  stable.  The  following 
example  illustrates  this. 


Example  1 


For  A  = 


1 
-2 


and  choosing  Q  =  I , 


it  follows  that 
C  =  -(AT  +  A) 


T 


The  system  is  stable,  but  C  is  indefinite.  Therefore,  x  x  is  not  a 
Lyapunov  function  for  the  system  described  by  the  given  transition 
matrix. 
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The  following  theorem  has  not  generally  appeared  in  the  literature 

and  is  formulated  here  as  a  new  result. 

Theorem  3.3      If,  for  a  stable  system,  its  transition 
matrix  A  is  real  and  commutes  with  its 
transpose  (i.e.,  AA  =  A  A)  then 

V  =  ||x||2  (i.e.,  Q  =  I) 

is  a  Lyapunov  function  for  the  system 
described  by  A. 

This  theorem  gives  a  sufficient  condition  for  the  square  of  the 

norm  of  the  states  to  be  a  Lyapunov  function.  To  show  this  we  take 

||x||  =  (x  x)  '    hence  V  =  x  x  and  show  that  A  +  A  must  be 

negative  definite.  But  A  +  A  is  a  real  symmetric  matrix  and  must 

have  real  negative  eigenvalues  [B2]  if  it  is  to  be  negative  definite. 

But  this  is  to  say  that  the  linear  system  described  by  the  transition 

matrix  A  +  A  is  stable.  Its  fundamental  matrix  is  then 

#  _  e(A  +A)t  ^  Nqw  if  ATA  =  AAT  then  we  can  write  ^     =  eA  t  eAt 

By  Theorems  2.1  and  2.2  of  Section  I  IB,  if   |$|  <_  k  then  the  system 

At 
is  stable.  But  since  A  is  stable   |e  |  <_  k,  so  we  write 

1 1$|  |  <_  k,k,  =  k 

For  discrete  systems  the  following  theorem  applies. 

Theorem  3. 3D     If  ADT  AQ  =  AD  ADT  then  a 

Lyapunov  function  for  the  discrete 
system  x(n+l)  =  AQ  x(n)  is 

V  =  xTx  with  aV  =  xT  (ADT  AQ  -  I)  x 

Suppose  there  exists  some  transformation  of  states  x  =  Sy 
such  that  the  transformed  system  matrix  D  =  S"  AS  is  normal 
(i.e.,  D  D  ■  DD  ) .  Then  we  can  apply  Theorem  3.3  and  a  Lyapunov 
function  for  the  transformed  system  is  V   =  y  y.  Since  stability 
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properties  for  linear  systems  are  invariant  under  a  similarity  trans- 
formation, we  may  apply  the  inverse  transformation  to  V  to  find  a 
suitable  Lyapunov  function  in  the  x  state  space.  The  transformation 
must  also  be  applied  to  V".  This  is  summarized  in  the  following  theorem, 

Theorem  3.4      If  a  non-singular  transformation 

exists  for  x  =  S"  y  such  that 
D  =  SAS"   is  normal,  then  a 
Lyapunov  function  for  the  system 
x  =  Ax   is: 

V  =  xTSTSx  with 

V"  =  xTST  (DT  +  D)  Sx 

=  xT  (ATSTS  +  STSA)  x 

A  similarity  transformation  may  be  applied  to  discrete  systems  as 

in  the  development  of  Theorem  3.4  for  continuous  systems. 


D.  LINEAR  FREE  SYSTEM 

Here  we  look  briefly  at  a  system  of  the  form  x  =  A(t)  x  . 
According  to  Theorem  3.1  of  this  chapter,  if  we  can  find  a  Lyapunov 
function  for  the  system  it  may  be  time  varying.  However,  we  may  also 
be  able  to  find  a  Lyapunov  function  which  is  not  explicitly  a  function 
of  time.  Assuming  that  we  can  find  a  positive  definite  0  such  that 
V  =  x  Qx  is  a  Lyapunov  function  then 

v  =  xT  [AT(t)Q  +  QA(t)]  x  =  -xT  C(t)  x 
In  order  to  apply  Theorem  3.1,  we  must  find  a  positive  function 
y( I | x | | )  such  that 

V(x,t)  <  -Y(| |x| I)  for  all  t. 

Apparently  we  must  somehow  remove  the  time  variation  of  A(t). 
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This  obviously  depends  upon  the  nature  of  A(t).  A  general  pursuit  of 
this  problem  is  not  very   fruitful.  However,  for  the  special  case 
considered  next  we  can  get  some  useful  results. 

E.  LINEAR  SYSTEM  WITH  PERTURBED  COEFFICIENTS 

The  special  case  of  the  linear  free  system  considered  here  is  that 
the  state  transition  matrix  A(t)  is  of  the  form  A(t)  =  A  +  G(t). 

T 

Now  we  assume  V(x)  =  x  Qx  and  compute 

tf(x,t)  =  -xTCx  +  xT(GT(t)Q  +  0G(t))x  (3.9) 

where  -C  =  ATQ  +  QA  . 

From  the  results  in  section  3C,  we  know  that  if  the  system  x  =  Ax 
is  stable  then  for  any  choice  of  a  positive  definite  matrix  C  we  can 
compute  a  positive  definite  matrix  Q.  Moreover,  if  G(t)  tends  to  zero 

as  t  increases,  then  there  exists  some  t  such  that  for  every     t  >   t 

o  J  o 

the  first  term  in  the  expression  for  V"(x,t)  will  dominate  and  v"  is 
negative  definite.  Hence,  the  assumed  quadratic  form  for  V(x)  is 
indeed  a  Lyapunov  function.  In  this  case  we  do  not  have  uniform  stabil- 
ity since  the  result  depends  upon  t  . 


Example  2 


Let   A(t)  = 


■1  + 


-2  + 


-1 

0 

Then 

A     = 

0 

-2 

G(t)  = 


Applying  Theorem  3C3  page  25  we  take  Q  =  I  then  C  = 


t 

0 

2   0" 
0   4 


and  C  will  dominate  for  all  t  ^  1/2.  Thus  this  system  is  ASIL  for 
any  tQ  >  1/2. 
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For  the  discrete  time  system 

x  (n  +  1)  =  Ax(n)  +  G(n)x(n) 
we  again  assume  a  Lyapunov  function  V  =  x  Qx.  Then  by  direct  compu- 


tation 


AV  =  V(x(n+1))  -  V(x(n))  =  -xT(n)Cx(n)  +  xT(n)QG(n)x(n) 


where  -C  =  ATQA  -  Q  . 

From  Theorem  3.2,  page  22,  we  know  that  for  any  positive  definite 
C  we  can  find  a  positive  definite  Q  if  A  is  stable.  Similar  comments 
apply  if  G(n)  tends  to  zero  as  for  the  continuous  case. 

F.   FORCED  LINEAR  SYSTEMS 

The  usual  application  of  Lyapunov' s  Second  Method  is  to  the  study 
of  the  dynamics  of  homogeneous  linear  systems  or  to  the  study  of  the 
dynamics  describing  the  motion  of  a  system  about  some  nominal  trajectory 
usually  attributed  to  the  forcing  function.  In  this  section  we  consider 
a  different  approach  which  leads  to  a  novel  interpretation  of  the 
Lyapunov  function  with  inputs  present. 

Consider  the  linear  system 

x  =  Ax(t)  +  Bu(t)  (3.10) 

Again  assume  a  Lyapunov  function  V  =  x  Qx  and  compute  its  time 

derivative. 

•  T     T  • 

V  =  x  Qx  +  x  Qx 

=  (uTBT  +  xTAT)Qx  +  xTQ(Ax  +  Bu) 
=  xT(ATQ  +  QA)x  +  uTBTQx  +  xTQBu 

V  =  -xTCx  +  2xTQBu  (3.11) 

From  Section  IIC  we  know  that  if  the  system  x  =  Ax  is  stable  we 
can  find  a  positive  definite  Q  for  any  given  positive  definite  C  such 
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that  -C  =  A  Q  +  QA.  Also,  note  that  the  bilinear  forms  uTBTQx  and 

x  QBu  are  transposes  of  each  other.  Since  they  are  scalars,  they  must 

be  equal,  and  they  may  be  combined  as  shown  in  the  second  term  of 

equation  (3.11). 

By  definition  2.1  in  Section  IIC  page  13  the  system  described  by 

(3.10)  is  stable  if  and  only  if  ||x||  <_  k, .  Thus  if  we  can  show  that 

V  given  by  (3.11)  is  negative  for  ||x||  >  k,  we  can  be  assured  that 

the  Lyapunov  function  will  converge  to  a  region  determined  by  ||x||  <  k. 

This  result  is  stated  in  the  following  theorem,  which  has  not  been 

stated  explicitly  in  the  literature  heretofore. 

Theorem  3.5      For  the  system  x  =  Ax  +  Bu(t), 

if  the  homogeneous  system  is  ASIL 
and  has  a  Lyapunov  function 
V  =  xTQx  with  V  =  -xTCx  and 
|  |u(t)|  |  <_  k2  for  all  t,  then  the 
system  is  stable  and  the  state  is 
sure  to  enter  a  region  defined  by 

2k2  MQII  ||B||  -  -§— 

||x||<.k1     =    }—        (3.12) 

Xl 

where     K%     ■     A  ,-„(C)     the  minimum 
1  mm 

eigenvalue  of  matrix  C     and     e  <  0 

with     |e|     arbitrarily  small 

Proof:     The  system  is  stable  by  definition  2.1   of  Section  IIC.     Take 

as  the  system  Lyapunov  function     V     =     x  Qx     then 

V     =     -xTCx  +  2xTQBu 

<  -xTCx+  2  ||x||  MQII  ||B||  ||u|| 
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also  since 


2    T 

A,   |x||    <_   X  Cx 


(3.12a) 


certainly  tf  <  -A-,  |  |x|  |  +  2k2|  |x|  |  ||Q||  ||B|| 

To  complete  the  proof  we  must  show  for  all  x  such  that   |x| |  .  k, 

then  V  <_  e  <  0 

When   Ixl  >  k,  ,  then 


x|  >_  k. 


x   > 


2k2||Q||  ||B 


and 
Thus 


A,||X 


>  2k2|  |x 


e  |  |  X 

K 


0  >  e  >_  A-i  1 1 x 1 1  +  2kp|  |x| 


B 


It  follows  that  v"  is  negative  when  |  |x|  |  >_  k,  ,  indicating  that  the 
systems  response  is  leaving  the  region. 

The  results  of  this  theorem  only  give  an  upper  bound  for  the  region 
to  which  the  state  will  converge  in  the  steady  state.  Unfortunately, 
this  bound  depends  upon  the  Lyapunov  function  chosen.  However,  the 
bound  is  determined  without  finding  the  solution  to  the  system  equations 


Example  3, 


Consider  the  stable  system 


""• 

— 

-      — 1 

-1 

1 

X 

+ 

1 

-1 

-1 

0 

sin  wt 


A  = 


-1 
-1 


B  = 


4                                                x  Cx 
Relation  (3.12a)  is  derived  using  the  Reyleigh  Quotient  

See  Bellman  [B3]  p.  110.  Ixl 
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A (A)  =  -1  ±  j 


AAT  =  ATA  = 


2    0 
0    2 


and  compute  a  steady  state  bound  for  the 
state  of  the  system. 

Theorem  3.3,  page  25  applies  to  the  homogeneous  dynamics,  therefore  a 

suitable  Lyapunov  function  is 


V  =  xTx 
that  is,    0=1 


\l     =  xT(A  +  AT)x  =  -xTCx 
~2    0 


C  = 


0 


and  the  homogeneous  system  is  ASIL.  The  results  of  Theorem  3.5  are 
used  to  compute  an  upper  bound  for   |x|   in  the  steady  state.  The 
parameters  for  equation  (3.12)  are  as  follows  taking  for  the  matrix 

||M||  =  [tr(MMT)]1/2     (see  Appendix  A) 

l|u(t)||   = 

linil   =  i 

Am.  (C)   =   2 
mm 


«=  (■„> 


sin  cot |  |  <  k2  =  1 


B 


=  1 


e  <  0  and  |e|  is  arbitrarily  small 


It  follows  that 


x|  |  </k,s 


2k2||Q||  ||B|| 


=  1 


It  should  be  noted  here  that  if  C  =  A, I  the  computed  bound  k, 
is  independent  of  x,.  This  follows  because  -C  =  A  Q  +  QA.  If  for 
C  =  I  the  resulting  Q  is  Q, ,  then  for  C  =  A, I  the  resulting  0  is 
A,Q,.  Thus   x,  cancels  out  in  the  expression  for  k, . 
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The  same  theorem  holds  for  discrete  systems  with  appropriate 
changes  in  notation. 

It  is  important  to  note  that,  for  linear  systems,  the  forcing 
functions  have  no  affect  upon  the  asymptotic  stability  properties  of 
the  system.  They  do  however,  influence  the  behavior  of  the  system 
in  the  steady  state. 

G.  TRANSIENT  ESTIMATION 

Since  the  level  curves  of  V  give,  in  some  sense,  a  measure  of 
distance  from  the  origin,  V  may  be  used  to  estimate  the  convergence 
rate  of  asymptotically  stable  systems.  The  following  is  taken  from 
Kalman  and  Bertram  [Kl]  with  some  minor  additions. 

Consider  the  obvious  inequality: 

V(x,t)  =  |||^|  V(x,t)  1-n-|  V(x,t)  (3.13) 

where  n-,  is  the  minimum  of  the  ratio  -  nrrrr*? )     in  some  region 
(radius  r)  of  state  space  excluding  the  origin.  In  cases  where  the 
basic  stability  theorem  applies  (Theorem  3.1  page  21)  we  may  take 

;   0  <  ||x||  <_  r\  (3.14) 

Integrating  equation  (3.13)observe  that 


inf 
x 


y! 

X 

) 

b( 

X 

) 

t 

V         ,1+ 

v     dt     <     - 

t 

t 

t 

0 

0 

nl 


dt 


9  ^(y)  >  0  *  llyll  1  rJ>  implies  the  greatest  lower  bound 
of  f(y)  for  all  y  such  that  0  <  | |y|   _  r. 
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V(x(t),t)  lV(xo,to)  e"nT(t-to)  (3.15) 

The  number,  -~  ,  is  an  upper  bound  for  the  time  constant 
describing  the  convergence  of  the  Lyapunov  function  following  the 
system  trajectory.  If  the  system  is  linear  the  Lyapunov  function  is 
a  quadratic  form  and  -^-      is  the  corresponding  bound  for  the  con- 
vergence of  the  system  trajectory. 


Similarly  we  can  write  V"  >_  -n?V  where  n?  1S  a  maximum  of  the 
o  — w-  ,  giving  a  lower  bound  for  the  system  time  constant. 
Kalman  [Kl]  has  carried  out  the  minimization  to  obtain  n- 


'1 
fo  the  linear  time  invariant  case.  These  results  are  stated  here. 

If  the  system  x  =  Ax  is  ASIL  and  has  a  Lyapunov  function  V  =  x  Qx 

with  8  =  -x  Cx  we  take 

n]     =        mxin(xTCx,  V=l)  (3.16) 

The  minimization  leads  to 

-Similarly 

n2     =     A^CQ"1)  (3.18) 

Figure  1   illustrates  the  relationship  of     V(x(t))     and  those  functions 
formed  by  the  n  bounds. 

These  convergence  estimates  depend  upon  the  Lyapunov  function 
chosen.     However,  there  are  certain  cases  where  the  estimates  are 
identical   to  those  obtained  directly  from  the  system  matrix  A.     The 
following  new  result  is  stated  as  a  theorem  the  proof  of  which  uses 
a  result  indicated  in  Bellman  [Bl]  involving  the  eigenvalues  of 
rational   functions  of  several  matrices,  as  follows: 
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V(x(t)) 


FIGURE  1 

Lyapunov  Function  and  convergence  bounds 
following  the  system  trajectory  x(t) 


^V(x(G,t)) 


FIGURE  2 

Effect  of  system  parameters  G  on  V 
G  minimizes  V 
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Let  A  and  B  be  commutative  matrices  and 
f(x,y)  be  a  rational  function.  Then  the 
characteristic  roots  of  f(A,B)  are 
f(a^»bj)a  where  a,-  and  b-j  are  the 
characteristic  roots  of  A  and  B  in  an 
ordered  sequence  which  depends  upon  the 
matrices  A  and  B. 

Theorem  3.6      If  Theorem  3.3  of  Section  II IC  page  25 

anplies  and  Q  is  chosen  as  the  identity 

matrix  then  -n*  =  2a  .   (A)  and  the 
I     mm 

convergence  estimate  using  Lyapunov 
functions  is  identical  with  that  of  the 
actual  system. 

T 

The  proof  is  immediate  since  if  0  =  I  then  C  =  -CA  +  A  ). 

T  T 

Taking  B  =  A   in  the  above  result  and  noting  that  A. (A)  =  A. (A  ), 

we  have  that  a. (A  +  AT)  =  2x.(A).  And  the  Amin(C)  =  2xmin(A). 

The  utility  of  this  method  is  obviously  not  in  estimating  the 
convergence  of  linear  time  invarient  systems  since  one  eigenvalue  prob- 
lem is  replaced  by  another.  Its  utility  in  this  thesis  is  in  application 
to  time  varying  systems  as  illustrated  in  the  following  example. 


Example  4. 


Let  A(t)  = 


o 


0 


-{*$) 


It  is  desired  to  estimate  the  convergence 

rate  for  the  system  described  by  the  given  A(t). 

Write  this  as 


A(t)  =  A  -  G(t)  = 


-1 

0 


0 


-2 


I    ° 


0 


By  applying  the  lemmas  developed  in  Section  IV  C  which  follows,  it  can 
be  shown  that  a  suitable  Lyapunov  function  is  (see  example  6  Section  IV  C) 


V  -  xTx 


-xT(A  +  AT)x  -  xT(G(t)  +  GT(t))x 


-Y 


■xTCx 


where  -C  =  A  +  A 


(3.19) 
T 
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For  this  example 
Q  =  I 


2    0 

0    4 


mi n  J y( | |x| 
X  (,   V(x). 


W«  i 


=  2 


Theorem  3.6  page  35  applies  also,  hence 

n,  =  -2x  .  (A)  =  2 
1      mirr  ' 

Suppose  that  V  is  a  function  of  a  matrix  of  parameters  G,  as  for 
example,  the  G(t)  in  equation  (3.19).  Now  if  for  some  G  =  G  , 
V  is  a  minimum  with  respect  to  G  then  the  ratio  -n-   is  also  a 
minimum,  hence,  ru  ,  is  maximum  (see  eq.  (3.13)).  It  follows  that  if 
n-,  is  maximum  for  G  =  G   then  the  system  involving  G   converges 
faster  than  systems  involving  other  G's.  Figure  (2)  illustrates  this 
concept  of  a  minimum  ratio  -n-  (G)  which  also  appears  in  Section  V 
in  connection  with  the  derivation  of  the  optimal  filter  gain.  In  this 
case  it  is  possible  to  choose  G(t)  so  as  to  maximize  the  system  con- 
vergence rate. 

Corresponding  results  to  the  foregoing  are  established  by  Kalman 
in  [K2]  for  discrete  systems.  The  main  difference  is  that  a  recursive 
relation  for  V(x(n))  is  established  as 

V(x(n+1))  <.e"2nTV(x(n)) 
where  T  is  the  sample  period  and  n  is  the  corresponding  constant  in 
the  continuous  case.  Again  the  minimization  leads  to 

e'2niT   =   Wof1) 

mm 


36 


H.   CONCLUSION 

General  results  are  difficult  to  obtain  using  the  foregoing  methods 
for  time  varying  linear  systems.  Often  the  results  give  sufficient 
conditions  only,  and  upper  bounds  which  depend  upon  the  chosen  Lyapunov 
function.  However,  there  is  an  advantage  in  that  these  results  are 
always  obtained  easily  without  solving  the  system  equation. 

The  approach  of  restricting  the  Lyapunov  function  to  be  explicitly 
time  invariant  was  considered  and  serves  as  an  introduction  for  a  unified 
approach  to  the  special  system  problem  considered  in  the  next  chapter. 
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IV.  THE  HOMOGENEOUS  FILTER  OBSERVER  EQUATION 

A.  INTRODUCTION 

The  remainder  of  this  thesis  is  concerned  with  the  study  of 
the  basic  recursive  averaging  filter-observer  using  the  foregoing 
concepts  of  stability  theory. 

The  filter  problem  is  usually  formulated  by  considering  a  signal 

model  driven  by  white  noise  and  whose  output  measurements  contain 

additive  white  noise.  Such  a  model  is  described  by 

x  =  Ax  +  w(t)  (4.1) 

z  =  Hx  +  v(t) 

where  A  is  an  nxn  state  transition  matrix,  H  is  an  mxn  measurement 

matrix,  w(t)  and  v(t)  are  n  and  m  vectors  respectively  of  uncorrelated 
white  noise. 

Filtering  is  accomplished  by  recursively  operating  on  the 

ovservations  z  to  estimate  the  values  of  x.  The  estimates  are 

denoted  x.  The  filter  dynamics  are  given  by 

x  =  Ax  +  G(t)(z-Hx)  (4.2) 

where  G(t)  is  the  nxm  filter  gain  matrix  the  choice  of  which  determines 

the  filter  characteristics.  If  the  error,  e  =  x-x,  is  formed  and 

the  expected  value  of  the  error  squared  (P  =  E[ee  ])  is  minimized 

to  determine  G(t)  then  (4.2)  becomes  the  optimal  integral  squared 

error  filter  derived  by  Kalman  and  Bucy  [K3].  If  G(t)  is  determined 

by  some  other  method,  the  filter  is  suboptimal  (for  example  see  S2). 


The  optimal  gain  is  given  by 
G(t)  =  P(t)HTR_1 

where  P(t)  satisfies  the  matrix  Ricatti  equation 


G(t)     P(t)HTR_1  (4.3) 
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P  =  AP(t)  +  P(t)AT  -P(t)HTR_1HP(t)  +  W  (4.4) 

and     R  =  E[vv  ],  W  =  E[ww  ] . 

Suppose  we  rewrite  equation  (4.2)  so  that  it  is  in  the  form  of 
equation  (4.1).  We  obtain 

x  =  (A-G(t)H)x  +  G(t)z  (4.5) 

This  equation  is  interpreted  as  a  system  having  a  state  transition 

matrix,  A-G(t)H,  with  the  input  vector  z.  Thus,  we  are  indeed 

taking  the  observations  z  and  filtering  them  to  provide  the  output 

state  vector  x. 

Let  us  also  find  the  equation  for  the  error  dynamics.  Define 
the  error,  e  =  x-x,  to  be  the  difference  between  the  actual  states 
and  the  estimated  values  out  of  the  filter. 

e  =  x-x  =  Ax  +  w(t)  -  Ax  -  G(t)(z-Hx)         (4.6) 
Noting  from  (4.1)  that  z  =  Hx  +  v(t)  we  write 

e  =  Ae  -  G(t)[Hx  +  v(t)  -  Hx]  +  w(t)  (4.7) 

and     e  =  (A-G(t)H)e  -  G(t)v(t)  +  w(t)  (4.8) 

Observe  that  the  homogeneous  dynamics  for  the  filter  equation 
(4.5)  and  the  error  equation  (4.8)  are  identical.  Both  are  described 
by  the  state  transition  matrix  A  -  G(t)H.  Consequently  the  results 
of  this  section  apply  to  both  the  output  of  the  filter  and  to  the 
actual  error. 

Now  suppose  that  G(t)  is  a  matrix  of  constants,  G  .  Then  the 
filter  transition  matrix  becomes  A  -  G  H  and  the  problem  reduces 
to  an  ordinary  linear  time  invariant  case.  Moreover,  if  G  =  I 
the  filter-observer  equation  reduces  to  that  form  derived  by 
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Luenberger  [LI]  for  observation  of  states. 

In  discrete  time  problems  the  signal  model  and  optimal  filter 

equations  becomes  respectively: 

x(n+l)  =  $(T)x(n)  +  w(n)  (4.9) 

z(n)  =  Hx(n)+v(n) 

x(n)  =  *(T)x(n-l)  +  G(n)[z(n)  -  H*(T)x(n-l)]       (4.10) 

G(n)  =  P(n/n-l)HT[HP(n/n-l)HT  +  R]"1  (4.11) 

P(n/n)  =  P(n/n-l)  -TG(n)HP(n/n-l ) 
P(n+l/n)  =  <fP(n/n)$  +  W 

where  the  notation  P(n+l/n)  is  the  estimated  covariance  of  error 

at  time  n+1  given  the  previous  n  corrected  estimates. 

For  discrete  systems,  we  find  equations 

x(n+l)  =  ($(T)  -  G(n)H*(T))x(n)  +  G(n)z(n)        (4.12) 
e(n+l)  =  ($(T)  -  G(n)H*(T))e(n)  -  G(n)v(n)  +  w(n)   (4.13) 
Again  the  unforced  dynamics  for  the  filter  states  and  the  estimation 
error  are  the  same. 

Bona  [B4]  has  considered  the  problem  of  choosing  the  eigenvalues 
of  the  matrix  (*(T)  -  G(n)H$(T))  to  give  a  desired  convergence 
property  (namely,  exact  estimate,  of  the  observation  in  a  finite 
number  of  sample  periods  for  noiseless  systems)  leading  to  the 
specification  of  G(n).  He  has  found  that  the  observer  can  be  made 
to  converge  as  rapidly  as  desired.  However,  the  observer  is  highly 
susceptible  to  error  due  to  the  presence  of  noise. 


Luenberger  [LI]  has  shown  by  using  linear  transformation 
theory  that  unobserved  output  states  of  a  free  linear 
system  may  be  observed  by  passing  the  observed  output 
states  through  another  linear  system  whose  transition 
matrix  is  A-H,  provided  A-H  is  stable. 
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The  filter  gain  matrix  G(t)  or  G(n)  has  a  twofold  purpose  in 
the  filter  observer  equation. 

1.  To  provide  observable  outpus  of  unovserved  states. 

2.  To  provide  some  programmed  weighting  of  the  actual 
observations  in  an  attempt  to  compensate  for 
errors  introduced  by  noise. 

Further-properties  -of  G  are  studied"! n  the  following  sections. 

B.   FILTER  STABILITY 

Stability  of  the  optimal  filter  has  been  demonstrated  in  general 

by  Kalman.  A  more  general  statement  of  the  following  theorem 

appears  in  Kalman  and  Bucy  [K3]. 

Theorem  1      (Kalman)  [K3].  For  a  linear  time  invar- 
iant signal  model  also  having  stationary 
statistics,  if  the  model  is: 
(i.)   completely  observable 
(ii.)  completely  controlable 
then  the  optimal  filter  is  asymptotically 
stable. 

It  is  evident  that  even  if  the  signal  model  equation  (4.1) 
section  IV.  A  is  unstable,  the  filter  may  still  be  asymptotically 
stable.  That  is,  the  matrix  (A-G(t)H)  in  equation  (4.5)  has 
its  own  stability  properties.  Thus,  it  is  not  necessary  that  the 
signal  model  be  stable  for  the  filter  to  be  stable. 

Deyst  and  Price  [Dl]  have  derived  conditions  for  asymptotic 
stability  of  the  discrete  optimal  filters  by  following  the  previous 
work  of  Kalman  for  the  continuous  optimal  filter.  Their  method 
uses  a  time  varying  Lyapunov  function  namely  choosing  Q  =  P"  (n) 
where  P(n)  is  the  discrete  covariance  of  errors.  Then, by  applying 
the  conditions  for  control abi 1 ity  and  observability,  they  found 
the  desired  bounding  functions  needed  to  apply  the  basic  Lyapunov 
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stability  theorem.  The  techniques  developed  in  this  paper,  namely 
requiring  the  Lyapunov  function  to  be  time  invariant,  generally 
apply  to  a  more  limited  class  of  signal  models  (i.e.,  asymptotically 
stable  ones).  However,  the  method  is  also  much  simpler  and  applies 
to  suboptimal  filters  as  well.  Such  a  study  leads  to  insight  into 
what  the  filter  gain  must  be  to  insure  a  stable  filter  and  may  in 
some  cases  be  used  to  actually  choose  the  G  matrix  requiring  that 
the  filter  be  stable. 

Generally,  the  technique  proposed  here  may  be  outlined  as 
follows: 

Assume  that  there  is  a  Lyapunov  function  of  the  form  V  =  e  Qe 
then  compute 

V  =  eT(ATQ  +  QA)e  -  eT(HTGT(t)Q  +  QG(t)H)e        (4.14) 
In  some  cases  equation  (4.14)  may  be  used  directly  to  determine 

the  matrix  Q.  Whereas,  in  other  cases, Q  may  be  determined  from 
-C  =  A  Q  +  QA  by  applying  Theorem  3.2  page  22.  If  A  is  asymptot- 
ically stable  then  V  of  (4.14)  becomes 

V  =  -eTCe  -  eT(HTGT(t)Q  +  QG(t)H)e  (4.15) 
If  equation  (4.14)  is  negative  definite  for  t>  t   then  the 

filter  is  ASIL  by  theorem  3.1  on  page  22.  Negative  definiteness  of 

V  may  be  shown  for  all  t  >  t  .  Alternately  it  is  possible  to  find  a 

bounding  function,  y(||x||),  as  required  by  theorem  3.1.  Obviously, 

if  the  second  quadratic  form  in  equation  (4.15)  is  at  least  positive 

semi-definite  for  all  t  >  t   then,  V  <  -e  Ce  and  the  conditions 

o        — 

of  theorem  3.1  are  met. 
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The  fact  that  the  filter  is  ASIL  means  that  considering  the 
filter  as  a  linear  system  free  of  input  forcing  functions,  the 
error  starting  at  any  arbitrary  initial  value  will  tend  to  zero 
as  t  increases. 

This  tells  us  nothing  about  the  actual  estimates  or  error. 
These  questions  must  be  considered  with  the  input  present  which 
due  to  its  stochastic  nature,  limits  the  asymptotic  value  of  the 
effective  error.  Such  consideration  appears  in  Chapter  V. 

For  discrete  systems,  taking  the  same  Lyapunov  function  as 
for  the  continuous  case  and  using  Theorem  (3. 2D)  page  24,  we  obtain: 
AV(e,n)  =  eT0T(I  -  G(n)H)TQ(I  -  G(n)H)$  -  Q]e     (4.16) 
AV(e,n)  =  -eTce  -  eT($TQG(n)H$  +  $THTGT(n)Q$)e 

+eT$THTGTQGH$e  (4.17) 

and     AV  =  -eTCe  -  2eT(<DTQG(n)H$)e  +  eT$THTGTQGH$e       (4.18) 
when    -C  =  $  Q$  -  Q 

The  following  serve  to  illustrate  the  kinds  of  filter  problems 
considered  in  this  thesis. 

1 .  Filter  Gain  Constant 

Case  1      We  are  given  G(t)  =  G  and  G  is  a  matrix 
of  known  constants.  For  this  case  we  may 
take  V  =  e  Qe  and  rewrite  equation  (4.14) 
as 

V  =  eT((A-GcH)TQ  +  Q(A-GcH))e         (4.18a) 
Since  A,H  and  G  are  all  known,  we  may 
apply  Theorem  (3.2)  page  22  directly  to 
determine  Q.  An  algorithm  for  numerical 
generation  of  Q  is  derived  in  appendix  B 
for  the  discrete  case  (Theorem  3. 2D). 
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Case  2      Choice  of  G  so  that  the  filter  is  stable. 
Now  use  equation  (4.15)  for  V  and  determine 
Q  based  upon  the  signal  model  by  applying 
Theorem  3.2  page  22.  This  requires  that  the 
signal  model  be  asymptotically  stable  if 
a  suitable  Q  is  to  be  found. 

The  second  quadratic  form  in  equation  (4.15)  is  now: 

eT(HTGcTQ  +  QGcH)e  =  eTDe 

D  =  HTG  TQ  +  QG  H 
c  v   x  c 

If  G  can  be  chosen  so  that  the  matrix  D  in  (4.19)  is 
c 

at  least  positive  semi-definite(p.s.d. ) ,  the  filter  is  sure 

to  be  stable  since  the  first  quadratic  form  in  equation 

(4.15)  is  already  negative  definite. 

Example  5 

Consider  again  the  system  described  in 
Example  3,  page  30.  Where 

"  l" 


(4.19) 


A  = 


Q  ■  I 


-1   1 
-1  -1 


B  = 


C  = 


Let  H  =  [1  0] 


In  case  2,  it  is  desired  to  choose  G  so 
that  the  matrix,  D,  in  (4.19)  is  at  least 
positive  semi-definite,  which  is  a  suffi- 
cient condition  to  insure  that  the  filter 

be  stable.  In  order  that  H  and  G   be 

c 

conformable,  G  must  be  of  the  form 
c 

V 
92 


44 


G  H  = 
c 


then 

I2 
and  from  (4.19) 

D  = 


^2 
0 


The  condition  that  D  be  p.s.d.  is  that  its 
principal  minors  be  non-negative  [ref.  B2]. 
That  is 


2gq  >_0 
2 

-g2  i 


0 


Obviously  g2  =  0  and  g,  >_  0  satisfy  the 
conditions.  Therefore,  any  filter  with 


[:] 


and  g,  >_  0 


will  be  ASIL. 

Now  suppose  we  are  given 
1 


Find  the  convergence  bound  n,  •  From  (4.18a) 
V  =  eT(MTq  +  QM)e 


where 


M  =  A-G  H  = 
c 


-2   1 
-1  -1 


We  may  compute  the  Lyapunov  function  Q  from 
theorem  3.2  for  any  p.d.  matrix  C  such  that 


-C  =  m'q  +  QM 
But  in  this  case 


1   =  M  M  and  we  may  use 
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theorem  3.3,  page  25.  Thus 


Q  ■  I 


C  =  -(M  +  M1) 


+4   0 
0  +2 


Th 


e  convergence  bound, n, ,  is  given  by  (3.17) 


nl 


=  WC^>  =  "min^     =  2 


2.  Special  Time  Varying  Gain 

Let  the  filter  gain  matrix  be  time  varying  such  that 
lim  G(t)  =  0  (4.20) 

the  null  matrix.  Here  we  have  a  problem  very  similar  to  that 
discussed  in  section  III  E,  and  those  comments  apply  here.  Also 


if  the  matrix 


(HTGT(t)Q  +  QG(t)H) 
in  equation  (4.15)  is  at  least  positive  semi-definite  for  all 


(4.21) 


t  >  t  ,  then  we  may  not  need  the  requirement  of  (4.20).  Sufficient 
conditions  for  the  matrix  (4.21)  to  be  positive  semi-definite  are 
developed  in  section  IV  C  which  follows. 
3.  The  Optimal  Filter 

Substituting  the  optimal  gain  G(t)  =  P(t)H  R~   into 
equation  (4.15)  we  obtain 

V  =  eT(ATQ  +  QA)e  -  eT(HTR_1HP(t)Q  +  QP(t)HTR_1H)e  (4.22) 
If  the  signal  model  is  asymptotically  stable,  Q  may  be  determined 
for  any  p.d.  matrix  C  as  before  and  (4.22)  becomes 

V  =  -eT(C  +  HTR_1HP(t)Q  +  QP(t)HTR_1H)e  (4.23) 
Now  the  condition  for  stability  from  (4.23)  is  that 

f(y,t)  ■  yT(C  +  HTR_1HP(t)Q  +  QP(t)HTK_1H)y  >  0     (4.24) 
for  all  t.  Moreover,  if 
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min  f(y,t)  >  0 
t 

then  the  condition  will  hold  for  all  t.  From  ordinary  calculus 

a  minimum  occurs  for  t  =  t-,  if 

df<y-V    >  o  and     d2f(y>V  >  0 

dt  ^r- 

Taking  the  derivative  of  equation  (4.24)  yields 

f(y,t)  =  yT(HTR"1Hf»(t)Q  +  Qp(t)HTR"1H)y  =  0      (4.25) 
One  way  for  (4.25)  to  be  zero  is  for 

P(t)  =  0  (4.26) 

Kalman  [K3]  has  shown  that  condition  (4.26)  yields  an  asymptotic 
or  steady  state  value  for  P(t)  in  the  limit  as  t  approaches  infinity 
by  solving  the  Ricatti  equation  (4.4)  page  39  for  P  with  P  =  0. 
However,  for  this  condition,  lim  P(t)  is  also  zero,  hence,  the 
foregoing  does  not  in  general  lead  to  the  desired  minimum.  Only 
in  the  trivial  case  where  P(t)  is  diagonal  does  this  method  lead 
to  the  desired  stability  condition,  since  the  diagonal  elements  of 
P(t)  are  monotonic  decreasing  to  the  limit. 

The  elements  of  the  matrix  in  (4.25)  are, in  general,  linear 
combinations  of  all  the  elements  of  P(t).  There  may  be  some 
t  =  t,  which  will  make  (4.25)  zero  other  than  in  the  limit  as 
t  approaches  infinity.  However,  t-,  is  generally  difficult  to 
find,  even  for  a  specific  example.  Often  it  is  necessary  to  compute 
the  trajectory  for  all  elements  of  P(t);  then,  using  this,  study 
the  properties  of  equation  (4.24)  at  various  values  of  time  in  an 
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effort  to  determine  t-, .  Besides  being  difficult,  the  method  ceases 
to  be  useful,  since,  once  P(t)  is  known  for  all  t,  complete  knowledge 
of  the  filter  is  at  hand. 

Further  research  is  required  in  the  area  of  finding  a 
suitable  bounding  function,  y(||e||),  for  equations  (4.22  or  4.23). 

For  suboptimal  filters  one  might  choose  G(t)  in  equation 
(4.15)  in  order  to  establish  the  stability  condition.  The  next 
section  attempts  to  build  some  insight  into  the  nature  of  G(t) 
required  to  insure  a  stable  filter. 

C.  STABILITY  CONSTRAINTS  ON  G(t) 

In  the  last  section  we  have  found  for  the  filter  error  dynamics 
that  if  V  =  e  Qe  is  a  Lyapunov  function  then 

V  =  -eTCe  -  eT(HTGT(t)Q  +  QG(t)H)e  (4.27) 

V  =  -eTCe  -  2eT(QG(t)h)e  (4.27a) 
where    -C  =  A  Q  +  QA  and  Q  may  be  determined  if  the  signal 

model    x  =  Ax  is  ASIL.  Equation  (4.27a)  is  written  as  shown 
for  convenience  even  though  the  second  term  is  not  a  true  quadratic 
form  (i.e.,  the  weighting  matrix  is  not  necessarily  symmetric). 
In  this  section  some  sufficient  conditions  are  found  which  insure 
that  V  in  (4.27)  be  negative  for  all  t.  If  these  sufficient  conditions 
hold  for  a  given  G(t)  then  the  filter  is  shown  to  be  ASIL.  On  the 
other  hand,  one  may  choose  G(t)  in  order  to  make  the  conditions  apply, 
thus  designing  a  filter  which  is  known  to  be  ASIL. 

In  applying  the  basic  stability  theorem  of  section  III  B,  we 
must  find  a  positive  function  y   such  that 
V(e,t)  :   -  Y(| |x| I)   for  all  t 
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Equivalently,  we  may  use  a  quadratic  form  for  the  bounding 
function  of  tf,  that  is,  we  must  find  a  positive  definite  matrix, 
D,  such  that: 

v"(e,t)  <  -  eTDe  =  -y  (4.28) 

Applying  this  to  equation  (4.27)  we  require: 

eT[C  +  HTGT(t)Q  +  QG(t)H]e  >  eTDe  (4.29) 

At  this  point  we  introduce  inequality  notation  for  symmetric 
matrices.  Namely,  A  >  B  means  that  (A-B)  is  positive  semi -definite. 
If  equality  is  removed,  then  (A-B)  is  positive  definite.  Writing 
relation  (4.29)  in  this  notation: 

C  +  HTGT(t)Q  +  QG(t)H  >_D   for  all  t  (4.30) 

The  matrix  D  is  important  if  one  is  to  use  the  technique 
described  in  section  3G  to  estimate  transient  response  of  the  filter. 
The  form  y  =     e  De  replaces  tf  in  the  ratio  vVV  used  to  determine 
n ,  hence  estimating  the  convergence  rate  of  the  filter.  Large  D 
means  that  the  gradient  of  y,  is  steep,  hence,  the  convergence  bound 
is  better  since  the  actual  filter  converges  faster  than  the  estimate 
given  by  n  •, .  Thus,  if  G(t)  is  known,  then  it  is  desired  to  choose 
the  largest  D  which  will  satisfy  (4.30).  On  the  other  hand,  if  we 
wish  to  design  the  filter  to  converge  faster  than  some  given  bound, 
then  D  is  given  and  (4.30)  provides  a  constraining  relation  for 
G(t). 

If  we  rewrite  (4.30)  and  apply  the  properties  of  norm  to  both 
sides,  we  get  another  bounding  relation  for  G(t). 

HTGT(t)f)  +  QG(t)H  >  D-C  (4.31) 
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,T,J 


|H'G'(t)Q| |  +  | |QG(t)H| |     D-C 


2||Q||  ||G(t)||  ||H||>  | | D-C i 


G(t) 


D-C 


01 


for  all  t, 


H 


(4.32) 


Intuitively  equality  in  (4.31)  and  (4.32)  should  result  in  a 
minimum  condition  of  the  left  hand  sides  for  all  t.  That  is  in 
(4.32) 


min  | |G(t)| | 
t 


D-C 


(4.32a) 


H 


A  similar  interpretation  of  the  equality  statement  for  (4.31) 
is  not  so  easily  conceived. 

The  basic  requirement  for  stability  given  by  (4.30)  or  (4.31) 
is  that 

min  Amin  (HTGT(t)Q  +  QG(t)H  -  D  +  C)  =  0  (4.33) 

Application  of  (4.33)  is  very  difficult.  However  sufficient  conditions 

to  insure  that  (4.33)  does  hold,  may  be  found  by  applying  an  important 

theorem  involving  bounds  on  the  eigenvalues  of  a  matrix.  The  Hadamard- 

Gerschgorin  theorem  is  stated  here  as  found  in  Bodewig  [B3]. 

Theorem  (Hadamard-Gerschgorin). 

The  eigenvalues  of  the  complex  matrix 
A  =  [a,-,-]  lie  inside  the  closed  domain 
G  consisting  of  all  circles  K. (i  =  l,...,n) 


with  centers  a. 
r. 


n  ui   ceri  iers>      d. 


and  radii     r.     where 


w 

If  we  apply  the  H-G  theorem  to  a  real   symmetric  matrix  and 
require  that    a,-..-  >/   |a..|     for  all   i   then  we  can  get  a  sufficient 


condition  test  for  positive  semi-definiteness. 
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Lemma  1 .   Let  A  be  a  real  symmetric  matrix  and  if 


iii  -  Ci 


*rj*|a..|  for  all  i,  then  A  is  positive 

semi -definite. 

The  proof  is  straight  forward  application  of  the  H-G  theorem. 

Since  A  is  real  symmetric  its  eigenvalues  (x(A))  must  be  real. 

Then  applying  the  H-G  theorem,  the  eigenvalues  of  A  must  lie  between 

min  (a..-  r . )  and  max  (a..  +  r. )  where  r.  =  /_.|a..|  as 
.     n   i       .    n    i         1     j?i  1J 

illustrated  in  Figure  3a.  But  the  hypothesis  requires 


implying 


a.  .  >  r.  =  *  i  fa 


a. .  -  r.  >_  0   for  all  i 


In  particular,  min  (a..  -  r , )  >_  0, 

which  implies  that  the  eigenvalues  of  A  must  be  non-negative.  Hence, 
A  must  be  positive  semi -definite  and  the  lemma  is  proved. 

Two  lemmas  will  now  be  proven  which  will  provide  a  method  by 
which  we  may  determine  the  best  bounding  matrix  D. 

Lemma  2.   Let: 

a.  M(t)  be  a  real  symmetric  matrix  of  elements 
m1;j  (t) 

b.  m..(t)  -  ZZjaUfCt)   =  a,(t) 

c.  a.(t)  >^  0   for  all  i  and  t 

Then  M(t)  is  positive  semi-definite  for  all  t. 
The  proof  is  immediate  from  Lemma  1.  Figure  3b  illustrates 
the  a.(t). 
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Im[x] 


All      \.     must  lie 


Re[x] 


Im[x] 


FIGURE     3a 

Hadamard-Gerschgorin  Theorem  applied  to 
a  real   symmetric  matrix 


A  =   (a..) 


n 


j^i 


l*,j 


^•(t) 


Re[\] 


FIGURE  3b 
Illustration  of  the  a-(t)  of  Lemma  2 
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As  a  result  of  lemma  2,  if  from  equation  (4.31)  we  set 

M(t)  =  HTGT(t)Q  +  QG(t)H  (4.34) 

then    M(t)  >_  D-C   for  all  t  (4.35) 

If  inf  a.(t)  =  0  for  any  i  then  then  condition  for  stability 
t   ' 
is  obtained  by  choosing  D  =  C.  However,  if  a.(t)  >  0   for  all  t, 

then  we  should  be  able  to  choose  some  D  >  C  giving  rise  to  a  tighter 

restriction  on  G(t) . 

Lemma  3.   In  Lemma  2  let: 

a.  a. (t)  =  min  a-  (t)   for  all  t 

i 

b.  a,  (t-, )  =  min  a.  (t)  >  0 

c.  m..(t)  -  mii(t])  >_  zZ^l^jCt)  -  mij(t1)| 

Then  M(t)  ^M(t-,)  with  equality  only  at  t  =  t-, . 

The  proof  of  Lemma  3  is  by  construction.  Choose  t,  as  specified 
in  conditions  (a)  and  (b),  thus  insuring  that  the  smallest  eigen- 
value of  M(t)  will  be  ppsitive  for  all  t.  For  the  matrix  (M(t)  -  M(t-,)), 
use  (b)  in  Lemma  2  to  compute  condition  (c)  in  Lemma  3,  and  the 
lemma  is  proved. 

It  is  apparent  that  the  conditions  in  Lemma  3  may  be  used  to 
insure  that  equation  (4.33)  holds.  That  is,  by  choosing  M(t)  as 
in  (4.34)  and  the  matrix  M(t-,)  as  in  Lemma  3,  we  may  rewrite 
equation  (4.33)  as  an  inequality  since  Lemma  3  only  specifies  a 

4 

boundary  for  the  eigenvalues  of  M(t).  The  result  is 

if      min  a  .  (M(t)  -  D  +  C))  >  0  (4.36) 

nnnv  v  '        — 

then     A  .  (M(tn)  -  D  +  C))  >  0  (4.36a) 

minv  v  1 '  ii~  K  ' 
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Now  if  Lermia  3  applies,  then  we  are  assured  that  (4.36a)  is  true 

for  some  choice  of  the  matrix  D.  However,  we  may  also  apply  Lemma  3 
directly  to  (4.31).  This  results  in 

M(t)  >     D  -  C  (4.37) 

Since  C  is  given  and  D  may  be  computed  directly,  by  replacing  the 
left  side  of  (4.37)  with  its  minimum,  M(t, )  ,  to  establish  the  equality 

M(t])  =  D  -  C  (4.38) 

This  follows,  since  by  Lemma  3,  M(t)  >_M(t, ). 

On  the  other  hand,  if  D  is  given,  as  established  by  some  con- 
vergence bound,  n-i  >  then  the  conditions  in  Lemma  3  may  be  applied  to 
(4.37),  thus,  establishing  constraining  relations  on  the  elements  of 
M(t). 

We  are  really  interested  in  the  matrix  G(t).  From  (4.34)  it  is 
obvious  that  every  element  of  M(t)  will  in  general  be  a  linear  combi- 
nation of  all  the  elements  of  G(t).  This  follows,  since  H  and  Q 
are  determined  by  the  signal  model  and  are  constant  matrices.  Thus, 
in  general 


M(t)  =  [mkl]  -  BCbu  9ij  (t))kl] 


(4.39) 


ij 


This  may  not  be  difficult  to  apply  to  specific  examples  since  H  and 
Q  are  often  very  simple  matrices  and  many  of  the  b.  .  =  0. 
The  following  illustration  refers  to  Example  4,  page  35. 
Example  6.       As  in  Example  4,  consider  the  system  described  by 

the  transition  matrix 

0 


A(t)  = 


1 


-0+-H 


o 


-(2+-J-) 


54 


This  is  equivalent  to  a  filtering  problem 
where  the  signal  model  is  described  as 
follows: 


The  system:   X  = 


r-1 

0 

r  "i 

W-j 

X 

+ 

0 

-2 

w2 

L.     - 

and  observation 


z  = 


where  A 


-1    0 
0   -2 


1 

0 

"vi" 

X 

+ 

0 

1 

v2 

*— 

H  =  I 


and  w, ,  w2,  v, ,  v~  are  white  Gaussian  noise 
inputs.  The  homogeneous  filter  dynamics  are 
given  in  (4.8)  as 


e  =   (A  -  HG(t))  e 
where 


A(t)  e 


G(t)  = 


1 


By  using  Lemma  2  we  will  derive  the  Lyapunov 
function  given  in  equation  (3.19)  of  Example  4. 

First,  assume  a  Lyapunov  function  of  the  form 

V  =  eTQe 
then  by  equation  (4.15) 

V  =  -eTCe  -  eT(HTGT(t)Q  +  QG(t)H)  e 

-e  De 

T  T    T 

where  Q  is  determined  from  -C  =  AD  +  OA.  But  in  this  case  AA  =  A  A 

and  we  may  use  Theorem  3.3,  page  25  to  find  0  and  C  .  The  result  is 
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Q  =  I 


C  =  -(A1  +  A) 


OH 


Now  in  order  to  apply  Lemma  2  we  find  from  equation  (4.34) 

2 


M(t)  =  G'(t)  +  G(t) 


This  follows  since  Q  =  H  =  I  . 

Then  by  Lemma  2  we  must  have 

2 

al  =  a2  =  T 


0 


for  all  t. 


inf 


Obviously,  a  >  0  for  any  t  >  t  >  0.  Since  T  o(t)  =  0  in  the 
limit  as  t  approaches  infinity,  the  condition  for  stability  is  ob- 
tained by  choosing  D  =  C  .  That  is,  -e  Ce  is  a  suitable  bounding 
function  for  V  .  Hence,  V  is  negative  definite  and  the  filter  is 
ASIL  by  Theorem  3.1,  page  21  . 

Lemma  3  applies  to  problems  where  the  a. (t)  in  Lemma  2  are  greater 
than  zero  for  all  t.  Consider  the  slight  extension  of  the  above  exam- 
ple by  letting 

1 


G(t)  = 


By  equation  (4.34) 


M(t) 


And  in  Lemma  3 


+  1 


0 


1 


+  2 


2  + 


4  + 


al  =  2  +  T 


a2  =  4  +  — 
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Take  t,  in  the  limit  as  t  approaches  infinity, 
(a)  and  (b)  in  Lemma  3: 
aft,)  =  2  >  0 


Then  from  conditions 


It  follows  that 


Mtt^  = 


and  condition  (c)  Lemma  3  reduces  trivially  to  Lemma  2.  Using  equation 
(4.38)  we  may  compute  the  matrix  D 

n(t^)  =  d  -  c 


D  =  M{t})   +  C  = 


4    0 
0    8_ 

By  using  the  results  of  Section  III  G  we  may  compute  a  convergence 
bound  for  the  filter  employing  the  given  G(t). 

The  problem  of  using  the  foregoing  Lemmas  to  obtain  constraints  on 
the  elements  of  G(t)  giving  rise  to  a  class  of  suboptimal  filters  is 
deferred  to  Chapter  VI. 

The  foregoing  Lemmas  apply  to  discrete  systems  simply  by  replacing 
t  with  n  .  Now  M(n)  is  specified  from  equation  (2D)  in  Section  4B. 
Namely: 


M(n)  =  $TQG(n)H$  +  $THTGT(n)Q$ 


$THTGT(n)QG(n)H$ 


(4.40) 


However,  (4.40)  may  be  simplified.  This  follows  since  $  is  non- 
singular  and  appears  as  a  congruent  transformation  on  the  inclosed 
matrix  and  does  not  affect  the  sign  definiteness  of  the  inclosed  matrix 
[ref.  B2].  So  in  effect  we  may  take 


M(n)  =  QG(n)H  +  HTGT(n)Q  -  HTGT(n)QG(n)H 


(4.41) 
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D.  CONCLUSIONS 

In  this  section  we  studied  the  homogeneous  dynamics  of  the  basic 
filter-observer  equation  (4.5)  and  its  error  dynamics  equation  (4.8), 
Since  these  homogeneous  dynamics  are  identical,  the  results  of  this 
section  apply  to  both. 

The  basic  approach  in  this  section  is  to  assume  a  time  invariant 
Lyapunov  function  and  compute  its  time  derivative,  V,  subsequently 
showing  that  V  is  negative  definite  for  a  given  filter  gain  matrix 
G(t).  This  was  found  difficult  to  do  in  general  terms.  Often  it 
is  possible  to  determine  the  quadratic  form  for  V  based  solely 
upon  the  signal  model.  This  requires  that  the  model  be  ASIL. 
It  is  desirable  to  choose  the  Lyapunov  function  based  only  upon  the 
signal  model  if  one  wishes  to  compare  the  convergence  of  two 
different  filters  designed  for  the  same  signal  model  by  using  the 
estimation  techniques  described  in  section  III  G.  Under  these 
conditions  the  Lyapunov  function  for  the  two  filters  is  the  same, 
but  the  associated  \l   is  modified  according  to  the  particular  filter 
(see  equation  (4.15)) . 

The  techniques  developed  in  this  chapter  cannot  be  readily 
applied  to  the  optimal  filter.  This  is  left  as  an  area  for  addi- 
tional research. 

Several  constraining  relations  between  the  bounding  y  function 
required  for  stability  and  characteristics  of  the  gain  matrix  were 
developed  in  section  IV  C.  In  general  these  relations  are  difficult 
to  apply,  but  in  some  simple  cases  they  lead  to  a  method  of  deter- 
mining G(t)  to  insure  a  stable  filter  which  converges  faster  than 
some  desired  bound. 
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V.  THE  FORCED  FILTER  EQUATION 

A.   INTRODUCTION 

In  this  chapter  the  forcing  functions  are  included  in  the 
filter  dynamics.  From  the  study  of  linear  forced  systems  in 
section  III  F,  recall  that  the  forcing  functions  do  not  affect 
the  asymptotic  stability  of  the  system  dynamics.  They  do,  however, 
determine  the  steady  state  properties  of  the  system.  Such  is  also 
the  case  for  the  filter-observer  system,  and  these  properties 
are  investigated  in  this  chapter. 

The  error  dynamics  for  the  filter  given  in  equation  (4.8) 
are  repeated  here 

e  =  (A-G(t)H)  e  -  G(t)  v(t)  +  w(t)  (5.1) 

As  in  chapter  IV,  assume  a  Lyapunov  function  of  the  form 

V  =  eTQe  (5.2) 

then  by  direct  computation  using  (5.1) 

•T     T  • 

V  =  e  Qe  +  e  Qe 

=  eT(ATQ  +  QA)  e  -  eT(HTGT(t)Q  +  QG(t)H)  e 
+  2wT(t)Qe  -  2vT(t)GT(t)0e 

V  =  -eTCe  -  eT(HTGT(t)Q  +  QG(t)H)  e 

+  2wT(t)Qe  -  2vT(t)GT(t)Qe  (5.3) 

where  Q  is  determined,  as  before,  for  any  given  p.d.  matrix  C  by 
the  signal  model  which  is  required  to  be  ASIL. 

It  should  be  noted  that  the  only  difference  between  equation 
(5.3)  and  equation  (4.15)  is  the  addition  of  two  terms  involving 
the  forcing  functions  w(t)  and  v(t).  However,  this  particular 
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application  is  different  because  the  forcing  functions  are  random 
vectors.  To  consider  this  aspect  we  use  the  expectation  operator 
and  rewrite  equation  (5.3).  Since  Q  is  completely  determined  for 
any  choice  of  positive  definite  C,  we  write  this  equation  in  terms 
of  Q.  Applying  the  matrix  identity  (5.4)  to  (5.3) 

xTMy  =  tr[MyxT]  =  tr[yxTM]  (5.4) 

results  in 

V  =  tr[Q(A-G(t)H)eeT  +  eeT(A-G(t)H)TQ 

+  2QewT(t)  -  2evT(t)GT(t)Q]  (5.5) 

In  order  to  remove  the  random  nature  from  this  problem  we  take 
the  expectation  of  V  and  V.  The  forcing  functions  w(t)  and  v(t) 
are  uncorrelated  white  Gaussian  random  processes.  Consequently, 
the  estimate  error  is  also  a  Gaussian  random  process.  We  define 
the  error  covariance  matrix  P  =  E[ee  ].  Thus  the  Lyapunov 
function  (5.2)  becomes: 

E[V]  =  trE[QeeT]  =  tr(QP)  (5.6) 

The  expectation  operation  is  done  over  an  ensemble  of  Lyapunov 
functions  each  of  which  is  the  result  of  one  experiment.  Taking 
the  expectation  and  applying  (5.4)  to  (5.5)  results  in 

E[V]  =  tr[Q(A-G(t)H)P  +  P(A-G(t)H)TQ] 

+  2QE[ewT(t)]  -  2E[evT(t)]GT(t)Q]  (5.7) 


'The  expectation  operator  is  defined  as  E[x]  =  /   xf(x)  dx 
where  f(x)  is  the  density  function  of  x  -°° 

[ref.  Papoulis,  PI] 


60 


To  finish  this  problem  we  must  first  determine  what  E[ew  ] 

and  E[ev  ]  are.  Equation  (5.1)  has  a  solution  of  the  form 

t  t 

e(t)  =  *(tnt  )e(t  )  +  /   $(t,x)w(x)dx  -  /   $(t,x)G(x)v(x)dT 

10     0.  . 

o  Lo        (5.8) 

But  w(t),  v(t)  and  e(t  )  are  uncorrelated  so  that 

E[wTw]  =  E[vTv]  =  E[wTe(to)]  =  E[vTe(to)]  =  0  (5.9) 
Also  define 

E[w(x)wT(t)]  =  W6  (t-x) 

E[v(x)vT(t)]  =  R6  (t-x)  (5.10) 

where  6(t-x)  =  ^  ^ 

Now  using  (5.8),  form  ew  . 

e(t)wT(t)  =  $(t,to)e(tQ)wT(t)  +  /   $(t,x)w(x)wT(t)dx 

to 

-  /   $(t,x)G(x)v(x)w'(t)dx  (5.11) 

lo 

Taking  the  expectation  of  (5.11),  interchanging  the  order  of 


integration  and  expectation  and  using  (5. 9), (5. 10) 

/ 
t 


E[e(t)wT(t)]  =  0  +  /t  $(t,x)W6(t-x)dx  -  0 


o 
=  1/2  [W]  (5.12) 

o 

This  follows  from  the  sifting  property  of  the  Delta  function  . 
Similarly  we  can  find: 

E[e(t)vT(t)]  =  -1/2  [G(t)R]  (5.13) 


The  Delta  function  sifting  property  applied  to  an  integration 
limit  at  the  location  of  the  Delta  function  is  used  here. 

b 
(i.e.,   /   f(x)6(x-b)dx  =  1/2  [f(b)] 
a 
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Substituting  (5.12)  and  (5.13)  into  (5.7)  yields 

+  G(t)RGT( 

(5.14) 


E[V]  =  tr[Q(A-G(t)H)P  +  P(A-G(t)H)TQ  +  QW  +  G(t)RGT(t)Q] 


Note  that  if  Q  =  I  in  (5.14)  the  result  becomes  the  error 
dynamics  developed  by  others.  Specifically  Athens  and  Tsi  [Al] 
arrived  at  (5.14)  in  their  derivation  of  the  optimal  filter  and 
Sims  and  Melsa  [S2]  also  arrived  at  this  result  for  specific  optimal 
filters. 

For  discrete  systems  a  difference  equation  for  E[V(n)]  analogous 
to  equation  (5.14)  is  easily  obtained  by  using  equation  (4.13) 
to  form 

E[eTQe]  =  tr[QP(n)]  (5.15) 

similarly  as  for  the  continuous  case.  The  result  for  P(n)  is 

P(n)  =  (I-G(n)H)P'(n)(I-G(n)H)T  +  G(n)RGT(n)       (5.16) 

P'(n)  =  $P(n-l')$T  +  W 
where  R  is  the  covariance  matrix  for  the  measurement  noise  and  W 
is  for  the  noise  driving  the  signal  model.  A  derivation  of  equation 
(5.16)  also  appears  in  Sorenson  [S3]. 

B.  DERIVATION  OF  THE  OPTIMUM  G(t) 

In  effect  the  stochastic  Lyapunov  function  (5.6)  for  the  filter- 
observer  is  a  linear  combination  of  all  the  covariance  of  errors 
between  the  estimates  (outputs  of  filter)  and  the  actual  signal 
being  estimated.  The  covariance  of  error  gives  a  measure  of  how 
well  the  filter  should  be  expected  to  perform.  As  such,  the  Lyapunov 
function  (5.6)  gives  a  scalar  measure  of  the  performance  of  the 
filter.  It  is  desirable  that  the  covariance  of  error  be  minimum 
for  all  t.  Thus  the  Lyapunov  function  (5.6)  may  be  used  as  a  cost 
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criterion  to  be  minimized.  The  optimum  filter  is  then  formed  by 
choosing  G(t)  such  that  E[v]  is  minimized  for  all  t.  If  (5.14) 
is  negative  for  all  t  >_  t  ,  then 

E[V]  <  0  (5.17) 

and  the  filter  is  ASIL  since  E[v]  >  0  for  all  t.  Moreover,  if 
E[V]  is  negative,  then  E[V]  is  assured  to  be  a  minimum  if  the 
magnitude  of  its  time  derivative  is  maximum  for  all  t.  The  optimum 
G(t)  can  be  derived  as  follows,  using  the  concept  of  a  gradient 
matrix.  The  following  formulas  for  the  partial  derivative  of  a 
scalar  function  (tr)  of  several  matrices  with  respect  to  one  of 
the  matrices  can  be  proven  (Athans  and  Tsi  [Al]). 

3  tr(ABT)  =  3  tr(BAT)   =  A  (5.18) 

3  B  3B 

3  tr(ACAT)  =  2AC  (5.19) 

3A 

where  C  is  symmetric. 

Taking  the  partial  derivitive  of  (5.14)  with  respect  to  G(t) 

using  (5.18)  and  (5.19)  yields 

3  E[V]   =  0  =  -  2PHTQ  +  2G(t)RQ 
3G 

and  solving  for 

G(t)  =  PHTR-1  (5.20) 

This  is  the  identical  result  derived  by  Kalman  and  others 
[Al ,  K3,  S3]  for  the  optimal  filter.  To  complete  the  derivation 
of  the  optimal  filter  we  must  show  that  (5.20)  results  in  a  minimum 
for  (5.14).  Note  the  second  partial  derivative  of  (5.14)  with 
respect  to  G,  by  applying  (5.18)  and  the  fact  that  R  and  Q  are 
symmetric. 
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rE[V]    =  2QR  (5.21) 

If      QR  =  RQ  (5.22) 

R  and  Q  commute,  then  the  right  hand  side  of  (5.21)  is  positive 

definite  since  Q  and  R  are  positive  definite.  This  follows  from 

the  result  indicated  in  Bellman  [Bl]  given  on  page  33  of  chapter  III. 

Now  if  (5.21)  is  positive  definite  then  (5.20)  does  result  in  a 

minimum  for  (5.14).  But  if  (5.17)  holds,  then  the  magnitude  of 

E[v]  is  maximum  (for  G  in  (5.20)),  hence  E[V]  in  (5.6)  is  a  minimum. 

Moreover,  the  ratio 

-E[V]  (5.23) 

E[V] 

is  also  maximum,  and  the  following  important  conclusion  is  reached. 

(l)  The  filter  using  the  6(t)  given 
by  (5.20),  is  assured  of  having 
the  most  rapid  convergence  rate 
to  the  minimum  covariance  of  error. 

This  follows  from  the  discussion  on  page  36  of  chapter  III,  section  G, 

The  foregoing  derivation  strictly  holds  only  if  conditions 

(5.17)  and  (5.22)  are  true.  However,  the  result  indicated  in 

(1)  above  does  give  engineering  insight  into  the  nature  of  optimal 

filters  and  to  desired  properties  of  suboptimal  filters.  That  is, 

it  is  possible  to  make  an  engineering  trade-off  between  convergence 

rate  and  steady  state  covariance  of  error  in  the  design  of  suboptimal 

filters  which  are  much  easier  to  implement  than  the  Kalman  filter. 

C.  CONCLUSION 

The  concept  of  a  Lyapunov  function  for  random  variables  was 
introduced  for  the  forced  filter  dynamics.  A  derivation  of  the 
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optimal  filter  gain  matrix  was  given  leading  to  the  conclusion  that 
the  optimal  filter  is  the  most  rapid  converging  to  the  minimum 
covariance  of  error.  As  a  consequence,  it  is  possible  in  the  design 
of  suboptimal  filters  to  make  an  engineering  trade-off  -between 
convergence  rate  and  steady  state  estimation  error. 
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VI.  AN  EXAMPLE  OF  A  SUBOPTIMAL  FILTER  . 

A.   INTRODUCTION 

Consider  the  discrete  form  of  a  system  described  in  phase  variable 
form  with  only  one  state  beinq  observed.  Namely: 


x  =  Ax 
y  =  Hx 


where 


A  = 


0  1  0 


0 


-a. 


0  1 
-a 


H  =  [1  0  ...  0] 
(6.1) 


It  is  assumed  that  A  is  stable  and  H  is  a  1  x  n  matrix,  and 
that  $(T)  has  been  found  as  described  in  Section  IIA  where  T  is 
the  sampling  period.  Note  the  transition  matrix  for  the  discrete 
filter  (4.12) 

x(n+l)  =  (*(T)  -  G(n)  H  *(T))x(n)  +  6(n)z(n) 

If  H  is  1  x  n,  then  G(n)  must  be  n  x  1  in  order  for  the  indicated 
multiplication  to  be  conformable.  This  is  indeed  fortunate,  since  G(n) 
is  chosen  to  insure  that  the  filter  is  stable  by  application  of  the 
lemmas  in  Section  IV  C,  and  these  lemmas  only  give  n  constraining 
relations  for  the  elements  of  G. 

Proceeding,  we  determine  the  system  Lyapunov  function  V  =  e  Qe 
for  aV  =  -e  Ce  by  finding  Q  such  that  $  Q$  -  0  =  -C  for  a  given 
positive  definite  C.  This  may  be  easily  done  using  the  algorithm 
developed  in  Appendix  B. 

Now,  in  order  to  insure  stability,  the  matrix  M(n)  in  (4.41)  must 
be  positive  semi -definite  in  accordance  with  Lemmas  1  and  2  developed  in 
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Section  IV  C.  The  argument,  (n)  ,  will  be  dropped  for  convenience 
in  the  following  and  the  notation  used  for  Q  is: 


Q  -  <►«>  n  x  n 

Repeating  equation  (4.41)  here 


M  =  QGH  +  HTGTQ  -  HTGTQGH 


and  for  G  =  (g.) 


j'  n  x  1 


(6.2) 


Multiplying  (6.2)  out,  noting  that: 
G  -  [g,  g2  .  .  .  gn]T 

GH=  [Gil'°]  nxn 


QGH  = 


£  b.  .  q . 
1   «  J  I 


HTST 


1 


0 


1     b  .  q. 


n  x  n 


I 


G  QG  is  in  fact  a  quadratic  form  in  the  vector  G  and 


HTH  = 


1  0  .  .  .  0 
0  0 


then  HTGTQGH  =  (  i     b.  .  g.  g.)  HTH 

Substituting  the  results,  (6.3)  through  (6.4),  into  (6.2) 
fl  b,.q.  -  E  b..q,q.)    I  b0  .q .  .  .  .  %   b  .g. 


M  = 


I  b2j9j 


E  b  .q. 


(6.3) 


(6.4) 


(6.5) 
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Now  the  a-  of  Lemma  2  are: 
1 


^  =  0 


n 

Z       b .  .  q  . 


i  =  2,3,  ...  ,n 


According  to  Lemma  2,  a.  >  0  however,  it  is  evident  that  for 


n  - 


i  >  2  we  must  choose  a.  =  0  since  they  can  never  be  positive.  This 
results  in  n  -  1  linear  equations  in  n  unknowns  depicted  by: 


n 

E 

J-l 


bij  9J  =  ° 


i  =  2,3, 


»n 


(6.6) 


Choosing  to  solve  these  n  -  1  equations  for  q2  ...  g   in  terms  of 
g,  ,  rewrite  equations  (6.6)  in  matrix  form 


b22   b23 


J32 


'n2 


2n 

—        - 
92 

93 

" 

nn 

% 

^m^l 

_         «. i 

-   b 


-  b 


21 


31 


-  b 


nl 


g1  (6.7) 


The  square  matrix  in  (6.7)  must  have  an  inverse,  since  it  is  a 

Q 

principal  minor  of  Q  which  is  positive  definite.    Thus,  (6.7)  may  in 
fact  be  solved  for  g2  ...  g   in  terms  of  g,.  The  remaining  relation 
necessary  to  determine  g,  is  given  by  a,  of  Lemma  2. 


l1 


J-l 


I,  .  u  . 

1j  yJ 


1SJ 


b.  .  g.  q.  -  0 
ij  yi  yJ 


(6.8) 


But  a.     must  be  non-negative,  giving: 


(6.9) 


g 

The  determinants  of  the  principal  minors  of  a  positive  definite 
matrix  are  positive  [ref.  B2]. 
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Using  equation  (6.7)  and  relation  (6.9)  a  family  of  G  matrices  may 

be  determined.  Since  some  a.  is  zero,  the  stability  condition  as 

discussed  in  Section  IV  C  following  Lemma  2  is  obtained  by  choosing 

D  =  C.  Referring  to  Section  III  G  the  convergence  bound  for  this  filter 

is: 

e"2nl   =  A  .  (DQ"1)  =  x  .JCn'1)  (6.10) 

rmnv  v  '  minv  •  '  v    ' 

But,  (6.10)  is  the  convergence  bound  for  the  signal  model  since  C 
was  chosen  for  the  signal  model.  This  means  that  the  transient  behavior 
of  the  filter  is  assured  to  be  no  worse  than  that  of  the  signal  model. 

Unfortunately,  the  foregoing  Lyapunov  technique  developed  cannot  be 
used  directly  for  optimal  filters  involving  the  phase  variable  form  of 
signal  model.  In  particular,  aV  ,  equation  (4.18),  has  not  been  shown 
to  be  negative  definite  for  the  optimal  filter  in  this  example.  How- 
ever, a  reasonably  good  comparison  of  the  suboptimal  filter  design  and 
the  optimal  filter  may  be  made  by  computing  tr(QP(n))  ,  equation  (5.15) 
for  both  filters  using  the  Q  determined  from  the  signal  model.  This 
follows  because  the  Lyapunov  function  (5.15)  is  a  scalar  function  of  all 
the  covariance  of  error  for  the  filter. 

B.  NUMERICAL  EXAMPLE 

In  order  to  summarize  the  foregoing  results,  a  simple  third  order 
numerical  example  is  given. 

A  suboptimal  filter-observer  is  to  be  designed  for  the  given  signal 
model.  The  filter  is  required  to  be  ASIL. 

The  signal  model  is  described  as  follows.  Given  the  signal  dynam- 
ics described  by  the  continuous  time  matrix  in  phase  variable  form 


69 


A  = 


0    1    0 

0   0    i 

-15  -17  -12 


(6.11) 


the  discrete  version  of  this  signal  model  for  a  sample  period  of  T  =  .15 
is  driven  by  a  white  Gaussian  random  sequence  whose  covariance  matrix  is 
given  as 


W  = 


0  0  0 
0  0  0 
0    0   0.0086 


(6.12) 


The  corresponding  discrete  system  for  (6.11)  was  found  for  T  =  .15 

by  using  a  digital  computer.  For  this,  a  common  NPS  subroutine  called 

AT 
"PHIDEL"  was  used  which  essentially  computes  the  series  for  e 

resulted  in: 


This 


$(.15)  = 


■ 


,994  .143  .0065 
0977  .884  .0653 
979     -1.207     .101 


(6.13) 


Only  one  state  of  the  signal  model  is  observed  and 

H  =  [1   0   0] 

The  observations  are  contaminated  by  additive  white  noise  whose 
covariance  matrix  is 


(6.14) 


R  = 


.5  0  0 
0  .5  0 
0     0     .5 


(6.15) 


This  completes  the  description  of  the  signal  model. 

The  first  step  in  the  design  of  a  suboptimal  filter  is  to  find  a 
Lyapunov  function  for  the  signal  model. 
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Assuming  that  C  =  I  the  Q  matrix  of  the  Lyapunov  function  for 
this  system  was  found  using  the  algorithm  developed  in  Appendix  B.  For 
this  problem  the  algorithm  converged  in  42  iterations  to  a  stopping 
criterion  of 


max 
ij 


d..(k) 


.0005 


(6.16) 


.th 


where  d.  .(k)  is  an  element  of  the  k    correction  matrix.  See 

'  J 

Appendix  B  for  further  explanation. 
The  resulting  Q  is: 


13.9 
Q  =         7.56 
0.357 

In  order  to  design  the  filter,  G  must  have  the  n  x  1  form 


7.56 

0.357 

13.6 

0.789 

0.789 

1.08 

(6.17) 


G  = 


q 


92 
93 


(6.18)  . 


and  equation   (6.7)  must  be  solved.     Substituting  from  (6.17)  into 


(6.7)  and  solving  for     g2     and     g-     results  in 

r  ~?   r    ' 

92 


13.6  0.789 

0.789         1.08 


93 


92 
93 


0.768 
-0.561 


-0.0561 
0.968 


-7.56 

-0.357 

P-7.56  " 
-0.357 


(6.19) 


g2     =     -0.56  g1 
g3     =       0.055  gj 


(6.20) 
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The  constraining  relation  (6.9)  is  used  to  obtain  a  restriction 
on  g,.  For  this  example  (6.9)  becomes 

13.9  g,  +  7.56  g2  +  0.357  g3  >  13.9  g^  +  2(7.56)  g,g2 

+  2(0.357)  gig3 
+  13.6  g22 
+  2(0.789)  g2g3 
+  1.08  g32 

Substituting  the  relations  (6.20)  and  reducing  this  equation  yields 

2 


9.69  g]  >  9.69  g^ 


1 


for  g1  > 


0 


(6.21) 


Thus  if    q,     is  constrained  by  equation  (6.21)  for  all     n   ,  and 
g2,g3     are  determined  by  (6.20),  the  filter  is  sure  to  be  stable  by 
its  construction  using  the  sufficient  conditions  of  the  Lyapunov  method 
in  conjunction  with  the  Hadamard-Gerschgorin  Theorem. 

The  simplest  kind  of  filter  to  implement  is  one  with  constant  gains, 


A  simple  time  varying  gain  which  satisfies   (6.21)  is: 

I  n 


,     for    n  =  1 ,2,3,. . . 


As  discussed  in  Section  VI  A,  a  comparison  of  these  filters  with 
the  optimal  filter  can  be  made  by  plotting  tr(QP(n))  for  each  filter. 

The  covariance  of  error  matrix  P(n)  was  computed  on  a  digital 
computer  for  the  optimal  filter  and  for  the  following  two  suboptimal 
gains  established  by  (6.20)  and  (6.21). 


Gl  = 


1 
-.56 
.055 


G2  = 


ri/n 
.56/n 
.055/n 


(6.22) 
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Computer  subroutines  for  these  calculations  appear  in  the  Computer 
Program  Section.  The  subroutine  "GAIN"  was  used  for  the  optimal  filter 
and  "SUBCOV"  for  the  filter  gains  in  (6.22).  Then  using  the  Lyapunov 
matrix,  Q  ,  in  (6.17) 

E[V]  =  tr(QP(n))  (6.23) 

was  plotted  on  the  same  axes  for  all  three  filters  as  shown  in  Figure  4. 

It  is  interesting  to  note  that  all  three  filters  are  stable,  and 
the  optimal  filter  does  converge  faster  than  the  suboptimal  filters  as 
predicted  in  Chapter  V,  Section  B.  The  optimal  filter  has  its  greatest 
advantage  in  the  first  one  or  two  sample  periods.  That  is,  the  correc- 
ted covariance  of  error  (P)  is  changing  more  rapidly  here  than  for  the 
suboptimal  filters. 

C.  CONCLUSIONS 

The  main  developments  in  this  thesis  show  that  it  is  feasible  to 
design  suboptimal  filters  which  do  perform  nearly  as  well  as  the  Kalman 
filter.  This  is  of  great  advantage  because  suboptimal  filters  often 
lead  to  a  much  simpler  implementation.  Such  filters  are  designed  so 
that  they  will  be  ASIL  and  converge  faster  than  a  known  bound.  The 
design  is  accomplished  by  using  Lyapunov  functions  and  requires  that  the 
signal  model  be  ASIL. 

There  are  two  areas  brought  out  by  this  thesis  which  require  further 
research. 

1.  Find  a  suitable  bounding  function  for  the  time  derivative  of 
the  Lyapunov  function  (4.22  or  4.23)  for  the  optimal  filters. 

2.  Find  bounding  relations  for  the  covariance  of  error  matrix,  P  , 
similar  to  those  developed  in  Theorem  3.5  for  the  system  state. 
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The  two  most  important  properties  of  the  filter-observer  system  are 
its  convergence  rate  and  its  ability  to  discriminate  against  noise.  By 
combining  any  results  obtained  for  1.  and  2.  above,  with  the  techniques 
developed  in  Chapter  IV  of  this  thesis,  a  design  procedure  for  sub- 
optimal  filters  could  be  formulated  which  would  insure  that  the  filter 
be  stable  and  also  give  information  about  its  two  most  important  prop- 
erties. Hence,  such  designs  could  be  formulated  and  compared  with 
other  such  designs  without  the  necessity  of  simulating  each  filter. 
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APPENDIX  A 
NORMS 

A  norm  is  a  scalar  function  of  several   variables  which  satisfy 

the  following  conditions.     The  norm  is  denoted     ||*||    . 

i)      |  |X|  |  ^  0     and     |  |X|  |   =0       iff     X  =  0 

ii)  ||uX||  =  | co |  ||X||       where  |u|  is  the  magnitude 

of  a  complex  scalar 

iii)  ||X  +  Y||<_||X||  +  ||Y||   Minkowski  inequality 
Possible  norms  for  an  n  vector  x 


|xtl         ||x||  =  ,???. 


1-1 


|x||  =  (  s  |Xi|2  )1/2  =  (tr(xxT))1/2 
|x||  =  (  I   \x/   )1/P       1  <  P  <  - 


Possible  norms  for  a  general  matrix  A 

||A||  =  [tr(AAT)]1/2         ||A||  =  i      |a..| 

l.j         J 

may         n  ?      ^2 

MAM      =       ax       t        la      I  ||All=         y        la      1^ 

|M|  |  .  l        | a.  .  |  |M|  h        | a. . | 

1   J-l    1J  ij    1J 

||A||  =  sup  ||Ax|| 
l|x||  =  1 

Other  properties  -  Schwarz  inequality 

I  I AB I I  <  I I At  I  x  I  I B| I       I I Ax I  I   <   I IAI I  x  I Ixl I 


Further  details  may  be  found  in  Ref.  Kl ,  SI,  Tl. 
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APPENDIX  B 
Lyapunov  Function  for  Linear  Discrete  System 

Davison  and  Man  [D3]  have  derived  an  algorithm  to  numerically 
generate  a  Lyapunov  function  for  linear  time  invariant  systems 
given  a  negative  definite  quadratic  form  for  V.  Here  a  similar 
algorithm  is  developed  for  use  with  linear  discrete  systems. 
In  section  III  C  it  is  shown  that  for  the  discrete  system: 

x(n+l)  =  ADx(n)  (1) 

a  positive  definite  Lyapunov  function  of  the  form: 

V  =  xTQx  (2) 

exists  with: 

AV  =  -xTCx  (3) 

being  negative  definite.  If  the  system  (1)  is  stable  Q  is 
determined  for  any  positive  definite  C  from  the  relation 

-C  =  ADQAD  -  Q  (4) 

In  general  it  is  difficult  to  solve  equation  (4)  for  Q  even 
given  that  C  =  I  let  alone  some  arbitrary  positive  definite  C. 
So  another  approach  is  taken  and  assumes  that  the  system  described 
by  (1)  is  ASIL. 

Proceeding,  we  write  the  Lyapunov  function  following  the 
system  trajectory  described  by  (1)  as: 

V(n)  =  V(0)  +  V(l)  -  V(0)  +  V(2)  -  V(l)  +  ...  +  V(n)  -  V(n-l) 

since    AV(k)  =  V(k)  -  V(k-l) 

n 
then    V(n)  =  V(0)  +   z   AV(i)  (5) 

i=0 
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Now  let  n  approach  infinity  then  V(n)  approaches  zero  by 

assumption  that  (1)  is  ASIL.  Then  (5)  becomes: 

V(0)  =  -    Z   AV(i) 
i=0 

Using  equation  (3) 

V(0)  =    z   x{  Cx.  (6) 

Noting  that  the  solution  to  (1)  is:  x.  =  Al  x    and  that 
V(0)  .  xJqxq 

We  rewrite  equation  (6): 

oo 

*°  QXo  =   j.j  Xo  (AD»TcAixo 
which  implies  that: 

oo      . 

Q  =    i   (An)  CA" 

i=0   u   u 

this  summation  will  converge  if  equation  (1)  is  stable,  since 

| x (AD) |  <  1   and  will  not  converge  if  (1)  is  unstable. 
The  proposed  algorithm  as  stated  here  is  easily  implemented 
on  any  computer  system  having  appropriate  matrix  subroutines. 

Dk  =  AD  Dk-1AD 

Qk+1  =  Qk  +  Dk         k  =  0,1,2,3,... 
with  initial  conditions: 

%  -  c  D-i  =  c 

and  termination  criterion: 

iiWi  <e 

The  matrix,  D,  is  the  correction  matrix  and  a  simple  conver- 
gence criterion  is  to  terminate  the  iteration  when  the  magnitude 
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of  all  elements  of  the  correction  matrix  are  less  than  some  desired 

amount.  That  is: 

max  |d..|  <_  e   where  D  =  (d..) 
ij  J 

Subroutine  "LYAP",  given  in  the  computer  program  section  following, 
was  used  to  generate  the  Lyapunov  Q  matrix  for  C  =  I  in  the  third 
order  example  given  in  chapter  VI.  The  result  is  given  in  equation 
(6.17). 
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THIS    SUBROUTINE    COMPUTES    THE    OPTIMUM    GAIN    MATRIX    AND 
THF    FRROR    COVARIANCE 


C 

c 


110 

111 

101 


108 


c 
c 


SUBROUTINE  GAIN(PK 
ODIMENSION  PKK(12,1 
1  ,G( 12,12),R(12,12 
2,TFMP1( 12,12),PHI( 
3TEMP2(12,12) 
G(K)  =  P(K/K-1)*HT 
CALL  TRANS(H,N,N,H 
CALL  PP0D(PKKM1,HT 
CALL  PROD(H,TEMPfN 
CALL  AOD(TEMPlf R.N 
CALL  RFCIP(Mf0.000 
IF  (KEP-2)  101,110 
WPITF(6,111) 
FORMAT  (5HKFR=2) 
CALL  PROD(TEMP,TFM 
NOTF  HERE  PKK(I  ,J 
o(K/K)  =  (I-G(K)*H 
CALL  PPOD(G,H,N,N, 
00  108  1=1, N 
00  108  J=l,N 
TFMP(  I,J)=-TEMP(I  , 
CALL  ADD  (HI, TEMP, 
CALL  PPOD(TEMP,PKK 
NOTE  HEPF  PKKMKI 
P(K/K-1)=  PHI*P(K- 
CALL  TRANSCPHI ,N,N 
CALL  PROD(PKK,PHlT 
CALL  PPOD(PHI,TEMP 
CALL  ADD(TEMPlfO,N 
RETURN 
END 


K,PKKM1,Q,R,PHI,H,N,M,G,HI,N0,M0,L0) 
2) ,0(12,12) ,H(12,12> 
) ,HI (12.12),HT(12,12),TEMP(  12,12), 
12,12)  ,PHIT(12,12),PKKMU12,12) 

*(H*P(K/K-1)*HT    ♦    R)**(-l) 

IT,ND,MD) 

",N,N,N,TEMPtND,MD,LD) 
i,N,N, TEMPI, ND,MD,LO) 
J.N. TEMPI, NO, MO) 
)001,TEMPlfTEMP2,KER,MO) 
1,101 

IP2,N,N,N,G,ND,MD,LD) 
I)     =    p(K/K)       WHERE 
)*P(K/K-1) 

N,TEMP,ND,MD,LD) 


J) 

N,N,TEMPfNO,MD) 
M1,N,N,N,PKK,ND,MD,LD) 
,J)    =    PCK/K-ll       WHERE 
1/K-1)*PHIT    ♦    Q 

,PHIT,ND,MD) 
,N,N,N,TEMP,NO,MD,LO) 
, N,N,N, TEMPI, ND, MO, LD) 
,N,PKKM1 ,ND,MD) 


SUBRO 
DIMEN 
£H( 12, 
#PHI( 1 
PK  =  ( I-G 
CALL 
CALL 
CALL 
CALL 
CALL 
CALL 
CALL 
CALL 
CALL 
PKM1=PHI 
CALL 
CALL 
CALL 
CALL 
RETUP 
END 


UTINE  SUBCOV 
SION  PKU2.1 
12).G(12.12) 
2,12),PHIT(1 

H)PKMi(I-GH) 

PROO(G,H,N,N 

SUB(HI ,TEM,N 

TRANSCTEM.N, 

PPOD(TEMfPKM 

PROD (D, PHI T, 

TRANS(G.N,N, 

PPOD(G,R,N,N 

PROD(TEM,D,N 

ADDCPK.PHIT, 

*PK*PHIT+W 

TRANS(PHI ,N, 

PRODCPK       ,PH 

PROD(PHI,TFM 

ADD(PKM1,W,N 

N 


(PKMlfW,R.PHI,H,N,M,G,HI,ND,MO,PK) 
2) ,PKM1 (12,12) ,W( 12.1 2 ),R( 12,12), 
.HI (12,12) ,TEM(12,12),D(12,12), 
2,12) 

T+GPGT 

,N,TEM,ND,ND,NO) 

,N,TEMfND,ND) 

N,PHIT,ND,ND) 

1 ,N,N,N,D,NO,ND,ND) 

N,N.N,PK,ND,ND,ND) 

D,ND,ND) 

,N,TEM,ND,ND,ND) 

,N,N.PHIT,Nn,ND,ND) 

N,N,PK       ,N,ND) 

^N.PHIT^D.NO) 
IIT,N,N,N,TEM,ND,ND,NO) 
SN,N,N,PKM1,ND,ND,ND) 
J,N,PKM1,ND,ND) 
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SURPOUTTNE    LYAP( SPHI , OL ,N, NO, MP, C) 

DIMFNSITN    SPHI (12,12)  ,QL (12,12 )  ,C( 12,12 ) , PH IT  (12,  12), 
«Ol( 12,12) .02(12,12) ,X2( 900) 
IY*PUN0V    FUNCTION 
OH    200A    1=1, N 
OH    200*    J=1,N 
0  2( I, J)=C( I, J) 

2004  0L( I, J)=C( T,J) 

CALL    TPANS(SPHI ,N, N, PHI T ,ND,NO) 

nn    200^    K=l,500 

CALL    PP0D(DHIT,D2,N,N,M,01 ,ND,ND,ND) 

CALL  PROP (01, SPHI ,N,N,N,D2,N0,NP,ND) 

CALL  APP(0L,P2,N,N,0L,N0,N0) 

TST=0. 

00  2006  1=1, N 

00  2006  J=1,N 

TTP=AB^(02( I  ,J) ) 

2006  CONTINUF 
IP(TIP.GT.TST)   TST=TIP 
IMTct. IE. 0.0005)  GO  TO  2007 

2005  CnNTINUF 

2007  WPITF(6,2008)  K 

2008  FOPMAT(«  ITTFRATION  K=»,I3) 
PPTUPN 

ENP 
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